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ABSTRACT 


This  report  summarizes  the  work  performed  during  the  first  year  of  a 
research  effort  to  determine  the  sound  fields  associated  with  jet  engine 
inlet  configurations,  A solution  approach  for  axisymmetric  bodies  based 
upon  the  integral  formulation  of  the  wave  equation  has  been  developed. 

This  solution  approach  circumvents  the  uniqueness  problems  which  normally 
occur  at  certain  frequencies  when  "straight  forward"  solutions  of  the  in- 
tegral equation  are  obtained.  A numerical  method  and  a computer  program 
for  solving  for  the  acoustic  field  associated  with  general  inlet  config- 
urations and  boundary  conditions  have  also  been  developed.  To  evaluate 
the  numerical  method,  computed  and  exact  results  are  compared  for  a sphere 
and  a finite  length  cylinder.  For  continuous  boundary  conditions,  the 
agreement  is  within  ten  per  cent  over  a range  of  nondimens ional  frequencies 
from  one  to  ten.  For  discontinuous  boundary  conditions,  the  numerical  errors 
increase  by  a factor  of  two.  This  report  presents  results  for  a given  inlet 
configuration  and  the  computed  and  exact  solutions  are  shown  to  agree  to 
within  ten  per  cent  over  the  nondimens ional  frequency  range  from  one  to 
ten. 
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I.  INTRODUCTION 


This  report  summarizes  the  results  obtained  during  the  first  year  of 
support  under  AFOSR  Contract  Number  F49620-77-C-0066.  This  contract  was 
initiated  on  February  1,  1977. 

The  research  conducted  under  this  contract  is  directed  towards  develop- 
ing analytical  techniques  for  predicting  the  characteristics  of  the  radiated 
sound  fields  from  jet  engine  inlets.  Such  capabilities  are  necessary  to  evalu- 
ate the  effectiveness  of  potential  sound  source  modifications  and  the  effi- 
ciency of  sound  suppression  techniques  for  fan  and  compressor  noise  attenua- 
tion in  inlets.  During  the  first  year,  the  conducted  research  efforts  have 
concentrated  on  the  development  of  an  efficient  analytical  technique  for  the 
prediction  of  the  radiated  fields  associated  with  lined  inlet  configurations. 
In  the  second  year,  experimental  investigations  will  be  conducted  to  provide 
data  for  comparison  with  the  theoretical  predictions. 

During  the  first  year,  a solution  approach  based  upon  an  integral  formula- 
tion of  the  wave  equation  has  been  developed  and  used  to  determine  the  charac- 
teristics of  the  sound  fields  of  several  previously  investigated  geometries. 
Efficient  numerical  techniques  have  been  devised  for  solving  the  integral 
equation,  and  the  necessary  computer  programs  have  been  written  and  tested. 
These  programs  are  now  capable  of  computing  the  surface  and  radiated  sound 
fields  for  arbitrary  geometries  with  lined  or  unlined  surfaces  and  sound 
sources  of  arbitrary  spatial  dependence.  These  capabilities  are  necessary  for 
the  investigations  of  sound  fields  from  jet  inlet  configurations. 

The  efforts  conducted  under  this  contract  has  resulted  in  three  publica- 
12  3 

tions  * ’ which  are  included  in  Appendices  A-C  of  this  report.  These  publi- 
cations provide  more  detailed  descriptions  of  the  research  efforts  conducted 
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under  this  contract.  The  research  performed  during  the  first  year  is  sum- 
marized in  the  following  sections. 

II.  ANALYTICAL  TECHNIQUE 

The  general  analytical  method  used  to  determine  the  radiated  sound  fields 
from  arbitrary  geometries  is  described  in  Appendix  A.  This  technique  is  based 
on  the  integral  form  of  the  solutions  to  the  wave  equation.  This  general  form- 
ulation has  been  specialized  to  axisymmetrlc  configurations,  which  are  applic- 
able to  jet  engine  configurations  which  are  of  interest  in  this  study,  in  Ap- 
pendix B, 

The  study  of  sound  radiation  involves  the  determination  of  the  acoustic 
field  over  an  infinite  domain.  However,  with  the  integral  formulation  of  the 
wave  equation,  the  acoustic  potential,  which  is  proportional  to  the  acoustic 
pressure,  can  be  computed  at  any  point  in  the  far  field  solely  from  the  values 
of  the  potential  distribution  at  the  surface.  Hius,  the  problem  is  reduced  to 
solving  for  the  acoustic  field  at  the  surface  only  instead  of  over  an  infinite 
domain. 

Several  problems  are  encountered  while  solving  the  integral  equation 
governing  the  surface  potential  distribution.  At  certain  frequencies  the  equa- 
tion fails  to  yield  a unique  solution.  These  frequencies  correspond  to  Internal 
eigenfrequencies  (or  resonant  frequencies)  of  the  geometry  under  consideration. 
This  nonuniqueness  manifests  itself  when  numerically  solving  the  integral 
equation  by  causing  the  coefficient  matrix  of  the  system  of  linear  algebraic 
equations  which  results  from  the  application  of  approximate  quadrature  to  the 
integral  equation  to  become  ill-conditioned,  causing  large  numerical  errors. 
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Using  a method  proposed  by  Burton  and  Miller  (Ref.  13  of  Appendix  A), 
this  behavior  can  be  eliminated.  This  method  consists  of  adding  the  inte- 


gral equation  for  the  normal  velocity  multiplied  by  a coupling  constant. 

It  is  then  proven  that  the  solution  for  the  acoustic  potential  field  from 
the  combined  equation  is  unique  for  imaginary  values  of  the  coupling  con- 
stant. This  analytical  method  was  therefore  incorporated  in  this  study. 

Although  other  techniques  can  be  used  to  avoid  the  uniqueness  problem  at  1 

certain  frequencies  (Refs.  4,  6,  and  14  of  Appendix  A),  the  combined  in- 
tegral equation  of  Burton  and  Miller  was  found  to  give  the  best  results, 
and  it  required  minimum  computation  times. 

In  order  to  use  the  Burton  and  Miller  method  two  problems  had  to  be 
resolved.  First,  a strong  singularity  exists  in  the  integrand  of  the  com- 
bined integral  equation  developed  by  Burton  and  Miller,  In  the  present 
study,  this  equation  was  reformulated  to  obtain  an  equation  containing  only 
weakly  singular  terms  which  could  be  handled  numerically.  The  second  prob- 
lem is  connected  with  the  choice  of  the  coupling  constant  used  by  Burton 
and  Miller  in  combining  the  integral  equations  for  the  potential  and  nor- 
mal velocity.  It  has  been  found  in  this  study  that  an  optimum  value  for 
this  parameter  for  use  in  numerical  computations  can  be  found.  Although 
Burton  and  Miller  showed  that  the  parameter  must  contain  a nonzero  imagi- 
nary component,  they  gave  no  indication  of  how  the  results  are  affected  by 
this  parameter.  The  value  which  gives  the  best  numerical  results  is  i/k 
where  i =\/^  and  k is  the  wave  number. 
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III.  NUMERICAL  METHOD 

A.  Integration  Procedure 

To  determine  the  acoustic  field  associated  with  a geometry,  the  integral 
equation  describing  the  surface  potential  distribution  must  first  be  solved. 
Using  this  distribution,  the  potential  at  any  exterior  point  can  then  be 
determined  to  generate  the  far  field  sound  pattern.  For  general  geometries, 
the  integral  equations  cannot  be  solved  exactly,  and  approximate  methods 
must  be  used.  These  methods  result  in  a system  of  linear,  algebraic  equations 
with  complex  coefficients  which  can  be  solved  by  complex  Gauss-Jordan  reduc- 
tion to  obtain  the  acoustic  potential  distribution  at  the  surface. 

For  the  axisymmetric  formulation  used  in  this  investigation,  the  surface 
shape  is  defined  by  a line  in  the  radial  r and  axial  z directions,  and  this 
line  is  rotated  about  the  axis.  The  integral  equation  can  then  be  separated 
into  two  line  integrals;  one  in  the  tangential  direction  and  one  along  the 
surface  contour  in  the  r - z plane.  In  the  tangential  direction,  the  line 
integrals  are  given  by  Eqns,  (15)-(17)  of  Appendix  B,  which,  in  general, 
must  be  solved  numerically,  A 96-point  Gaussian  quadrature  formula  was  used 
to  evaluate  these  integrals.  The  computational  error  is  approximately  in- 
versely proportional  to  the  number  of  points  used  to  evaluate  the  integrals 
in  the  tangential  direction.  Along  the  surface  contour  in  the  r - z plane, 
the  integration  of  Eq.  (18)  in  Appendix  B is  required.  The  integral  over 
the  perimeter  is  first  separated  into  integrals  over  n subintervals  of  either 
constant  or  varying  lengths.  The  acoustic  potential  is  assumed  constant  over 
each  subinterval  and  it  is  taken  outside  the  integral.  Finally,  a two-point 
Gaussian  quadrature  formula  is  used  to  evaluate  the  components  of  the  func- 
tions defined  by  Eqs,  (15)-(17)  of  Appendix  B in  the  r - z plane  over  each 
subinterval.  Increasing  the  order  of  the  Gaussian  quadrature  does  not  sig- 
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nificantly  affect  the  accuracy  of  the  computations.  However,  the  error  was 
found  to  decrease  proportionately  as  the  number  subintervals  n was  increased. 

In  another  solution  approach,  the  potential  was  assumed  to  vary  linearly 
over  each  interval  in  order  to  improve  the  accuracy  of  the  correct  potential 
values.  However,  the  results  using  this  linear  interpolation  scheme  were  not 
as  accurate  as  in  the  above-mentioned  approach.  The  cause  of  the  inaccuracies 

in  this  scheme  have  not  been  extensively  studied  because  of  time  restrictions,  | 

However,  the  errors  appear  to  arise  from  the  implementation  procedure  used. 

When  linear  interpolation  was  applied  to  a finite  cylinder,  problems  arose 
at  the  corner  points.  At  these  points  the  normal  to  the  surface  appearing 
in  the  integral  equation  is  undefined.  If  the  potential  is  assumed  constant 
over  each  subinterval,  subintervals  can  be  taken  on  either  side  of  this 
point  which  in  effect  avoids  the  corner  points.  The  method  used  in  applying 
the  linear  interpolation  about  these  points  strongly  influences  the  computed 
results.  Although  several  techniques  were  tried,  none  proved  entirely  satis- 
factory, Also,  for  general  surfaces  for  which  the  subintervals  may  be  of  un- 
equal length,  the  difficulty  in  implementing  the  linear  interpolation  tech- 
nique and  its  questionable  value  make  this  method  impractical.  Therefore,  it 
will  not  be  used  in  future  studies  unless  significant  improveiii'-  itj  can  be 
made. 

In  another  study  involving  the  numerical  evaluation  of  the  integral 
equations,  the  effect  of  the  coupling  constant  (a  appearing  in  Eq,  (23)  of 
Appendix  A)  was  investigated.  In  the  method  of  Burton  and  Miller,  the  inte- 
gral equations  for  the  acoustic  potential  and  normal  velocity  at  the  surface 
are  combined  into  one  equation.  The  terms  from  the  potential  equation  are  of 
order  k whereas  the  terms  from  the  expression  for  the  normal  velocity  are  of 
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2 2 
order  k . Thus,  as  the  frequency  is  increased,  the  terms  of  order  k 

dominate.  The  results  become  less  accurate  because  the  combined  equation 

in  effect  becomes  the  equation  for  the  normal  velocity.  This  equation, 

like  the  integral  formula  for  the  acoustic  potential,  yields  large  errors 

at  certain  frequencies  when  numerically  evaluated.  By  choosing  the  coupling 

2 

constant  to  be  i/k,  the  terms  of  order  k are  now  reduced  to  order  k.  Now, 
as  the  frequency  is  increased,  the  terms  from  the  expression  for  the  normal 
velocity  do  not  become  dominant,  and  the  uniqueness  problem  is  avoided  at 
all  frequencies, 

B.  Evaluation  of  Geometric  Parameters  and  Boundary  Conditions 

The  geometric  parameters  appearing  in  the  integral  Helmholtz  equation 
(see  Equation  (23)  of  Appendix  A)  are  the  distances  between  points  on  the 
surface,  the  normal  vector  at  each  point,  and  the  lengths  of  the  subinter- 
vals in  the  axial  plane.  The  method  for  computing  these  parameters  is  pre- 
sented in  Section  III  of  Appendix  C, 

There  are  two  types  of  boundary  conditions  which  must  be  specified  over 
the  surface.  The  first  consists  of  a forcing  function  which  generates  the 
acoustic  field.  In  a jet  engine  inlet,  most  of  the  acoustic  field  is  pro- 
duced by  disturbances  caused  by  the  interaction  between  the  flow  field 
produced  by  the  fan  blades,  and  the  stator  waves,  A literature  search  was 
conducted  to  determine  the  spatial  dependence  of  the  sound  generated  by  the 

stator-blade  interaction  so  that  the  resulting  radiated  sound  pattern  could 

4 

be  computed.  Because  of  the  complexity  of  the  resulting  expressions  , there 
was  not  sufficient  time  to  use  these  predictions  in  the  present  research  effort. 
However,  both  the  analytical  and  numerical  methods  used  in  computing  the  radi- 
ated sound  field  in  the  present  study  are  capable  of  handling  forcing  functions 
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of  arbitrary  spatial  dependence  in  the  r - z plane,  such  as  those  encountered 
in  jet  engine  inlets. 

The  second  type  of  boundary  condition  is  given  by  specifying  the  reaction 
of  the  surface  to  the  wave  motion.  For  rigid  surfaces,  the  normal  velocity 
(i.e„,  the  normal  derivative  of  the  acoustic  potential)  is  zero  and  all  the 
sound  incident  on  the  surface  is  reflected.  For  nonrigid  or  sound  absorbing 
surfaces,  the  normal  velocity  is  nonzero  since  the  surface  now  vibrates  in 
response  to  the  wave  motion.  The  normal  velocity  at  the  surface  is  propor- 
tional to  the  pressure  oscillations  (i.e.,  the  acoustic  potential)  of  the 
surrounding  fluid  and  the  constant  of  proportionality  is  called  the  surface 
admittance.  The  integral  wave  equation  involves  both  the  acoustic  potential 
and  the  normal  acoustic  velocity  which  means  there  is  one  equation  for  two 
unknowns.  By  using  the  admittance  relationship,  the  normal  velocity  can  be 
expressed  in  terms  of  the  acoustic  potential  and  the  admittance.  The  result- 
ing equation  can  then  be  solved  for  the  potential, 

Tlie  admittance  is  a measure  of  the  sound  absorption  characteristics  of 
the  surface.  In  jet  engine  inlets  the  surfaces  are  often  lined  with  Helm- 
holtz resonator  arrays  which  absorb  sound  and  reduce  the  noise  radiated  to 
the  surroundings.  Expressions  for  the  admittances  of  these  sound  absorbing 
devices  have  been  derived^  and  they  can  be  used  in  the  present  investigation. 
In  fact,  the  capability  exists  to  predict  the  sound  field  produced  from  a 
jet  engine  inlet  for  arbitrary  sound  source  and  admittance  characteristics. 

C.  Computer  Program 

A computer  program  written  in  Extended  FORTRAN  IV  has  been  developed 
for  use  on  a CDC  CYBER  70  computer  for  solving  the  system  of  linear  alge- 
braic equations  which  result  from  the  numerical  approximation  to  the  inte- 
gral wave  equation.  This  program  has  been  thoroughly  checked  out  using 
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simple  geometries  for  which  exact  solutions  can  b obtained.  The  program 
employs  standard  functions  common  to  all  FORTRAN  compilers  so  that  it  can 
be  used  with  minimum  modifications  on  other  computers.  Generality  is  main- 
tained in  order  to  accommodate  arbitrary  surface  geometries  and  boundary 
conditions.  In  the  cases  run  to  date,  the  computation  time  for  determining 
the  surface  potential  is  given  by  the  following  approximate  formula 

t = 0,05  (n)^ 

where  t is  the  computation  time  in  seconds  and  n is  the  number  of  subinter- 
vals used  in  the  numerical  evaluation  of  the  integral  equation.  Approximately 
one  second  per  point  is  required  for  the  far  field  potential.  For  the  cases 
run  thus  far,  the  run  times  have  been  from  20  to  140  seconds  for  20  to  53 
subintervals.  Efforts  toward  maximizing  the  programming  and  numerical  effi- 
ciency have  resulted  in  these  relatively  short  run  times. 


IV.  RESULTS  AND  SUMMARY 

A,  Simple  Geometries 

To  check  the  numerical  schemes  used  in  this  investigation,  preliminary 
computations  using  a sphere  and  a cylinder  of  finite  length  were  obtained. 
The  results  are  presented  in  Section  III  of  Appendix  B,  In  all  cases,  20 
subintervals  were  taken  in  the  r - z plane  and  a 20-point  Gaussian  quadrature 
was  used  to  evaluate  the  integrals  in  the  tangential  direction.  The  results 
for  these  simple  geometries  can  be  summarized  as  follows: 

(1)  The  coupling  parameter  used  in  the  Burton  and  Miller  (a  in  Eq,  (23) 
Appendix  A)  should  be  taken  as  i/k  where  k is  the  wave  number  and 
i is  \/^. 
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(2)  For  the  cylinder,  the  error  in  the  computed  results  increases  with 
increasing  frequency.  For  continuous  boundary  conditions  the  error 
is  less  than  107.  at  all  frequencies, 

(3)  Discontinuous  boundary  conditions,  where  the  admittance  is  specified 
over  part  of  the  surface  and  a forcing  function  over  the  remainder, 
decrease  accuracy  of  the  results.  The  computed  and  exact  values  agree 
to  within  107o  for  low  nondimens  ional  wave  numbers  (i.e,,  ka<  5 where 
a is  the  radius  of  the  cylinder).  Errors  of  407.  at  the  point  of  dis- 
continuity occur  at  a nondimens ional  wave  number  of  10,  The  remainder 
of  the  points  agree  to  within  127.  at  this  frequency, 

(4)  In  the  tangential  plane,  the  spatial  distribution  of  the  acoustic 
potential  varies  as  cos  m 0 where  m is  an  integer.  Increasing  m 
does  not  affect  the  accuracy  of  the  results  significantly, 

(5)  The  computed  far  field  acoustic  potentials  are  at  least  as  accurate 
as  the  computed  surface  potentials, 

(6)  The  far  field  results  are  accurate  at  distances  greater  than  the 
length  of  one  subinterval  from  the  surface. 

B.  Inlet  Configuration 

The  studies  of  the  acoustic  fields  of  the  sphere  and  cylinder  served  to 
evaluate  and  refine  the  numerical  procedures  and  programming  techniques.  The 
next  configuration  investigated  was  an  inlet  used  in  a study  by  NASA^.  This 
inlet  is  shown  in  Fig.  1 and  it  was  chosen  because: 

(1)  unlike  most  inlets  used  in  research  studies,  it  does  not  have  a bell- 
mouth  shape  but  is  shaped  like  a typical  inlet  used  in  existing 
aircraft; 
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(2)  complete  details  on  generating  the  inlet  boundary  are  given;  and 

(3)  it  is  being  used  in  a related  study  being  conducted  at  Georgia 
Tech  concerning  the  prediction  of  the  sound  field  inside  the 
inlet;  so  the  sound  field,  at  least  inside  the  duct,  can  be  com- 
pared with  results  obtained  independently  by  other  numerical 
methods , 

The  back  side  of  the  inlet  is  presently  assvimed  to  be  spherical. 

To  obtain  exact  results  for  comparison  with  the  numerical  computations, 
a spherical  source  was  assumed  to  be  placed  at  (r,z)  = (0,0).  The  acoustic 
potential  and  normal  velocity  for  this  source  can  be  readily  computed  at 
every  point.  In  particular,  they  can  be  computed  on  the  surface  of  the 
inlet.  The  value  of  the  normal  velocity  at  each  point  along  the  surface  of 
the  inlet  is  then  used  as  the  boundary  condition  in  the  integral  equation. 

From  this  boundary  condition  the  value  of  the  potential  can  then  be  numerical- 
ly computed  using  the  techniques  described  in  Chapter  III  and  compared  with 
the  exact  potential  known  from  the  spherical  source  solution.  As  seen  in 
Fig.  2,  the  normal  acoustic  velocity  distribution,  which  represents  a forcing 
function  is  highly  discontinuous  and  it  provides  a severe  test  of  the  numeri- 
cal techniques  employed. 

The  numerical  and  exact  solutions  for  the  surface  acoustic  potential  are 
compared  in  Fig,  2 for  32  and  54  subintervals  taken  along  the  perimeter  of 
the  inlet  in  the  r - z plane.  Because  of  the  errors  in  approximating  the  lengths 
of  each  subinterval,  the  exact  solutions  differ  slightly  as  the  distance  along 
the  perimeter  S increases.  The  centerbody  in  Fig.  1 extends  from  0<S<0.8, 
the  fan  inlet  covers  0,8  <S<  1,4,  the  interior  contour  extends  from  1 ,4  < S < 3,5, 
the  exterior  from  3,5  <S<  5,5,  and  the  circular  arc  lies  within  the  interval 
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Fi  ^,re  2.  Effect  of  Increasing  the  Number  of  Subintervals  in  Computing 
the  Surface  Potential  for  the  Inlet  Configuration  at  k a = 1 


5.5 < S<  7.45,  Increasing  the  number  of  points  decreases  the  error  pro- 
portionately as  indicated  by  the  data  in  Fig,  2 at  a nondimens ional  fre- 
quency to  a of  unity,  where  a is  the  radius  of  the  inlet  at  the  fan  entrance 
section.  The  absolute  average  error  in  the  results  decrease  from  10.2  per 
cent  for  32  subintervals  to  4.16  per  cent  for  53  subintervals.  The  computa- 
tion time  increased  from  53  seconds  to  143  seconds,  respectively. 

As  shown  in  Fig,  3,  the  errors  increase  with  increasing  frequency.  Like 
the  cylinder,  the  maximum  error  of  the  potential  for  the  inlet  configuration 
occurs  at  the  points  of  discontinuity.  The  average  error  increases  from  4,16 
per  cent  at  k a = 1 to  15  per  cent  at  k a = 10, 

For  the  data  in  Figs,  2 and  3,  the  acoustic  potential  is  assumed  constant 
in  the  tangential  plane.  The  results  for  a cos (m  0)  distribution  are  presented 
in  Fig,  4 at  ka  = 2.  These  results  show  the  insensitivity  of  the  accuracy 
of  the  computed  results  to  the  tangential  distribution  for  m = l,  2. 

Based  on  the  results  obtained  thus  far,  the  numerical  and  programming 
techniques  are  capable  of  yielding  reliable  results  for  arbitrary  geometries 
and  boundary  conditions.  At  higher  frequencies , (ka < 5)  it  appears  that  more 
points  must  be  taken  to  increase  the  accuracy  of  the  computed  results. 

Next  year,  experiments  will  be  conducted  to  measure  the  acoustic  field 
radiated  from  an  open-ended  pipe  for  comparison  with  the  computed  results, 

A parametric  study  of  the  effect  of  the  placement  and  quality  of  sound  treat- 
ment on  sound  abatement  in  an  inlet  configuration  will  be  conducted.  Further 
improvements  in  the  programming  and  numerical  methods  will  also  be  investi- 
gated. 
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Abstract 

This  paper  is  concerned  with  the  development  of  a procedure  for  generating 
the  sound  fields  radiated  by  arbitrarily  shaped,  three  dimensional  bodies  from 
an  integral  representation  of  the  solutions  of  the  Helmholtz  equation.  The 
method  of  Burton  and  Miller  is  employed  to  eliminate  the  nonuniqueness  in  the 
external  Helmholtz  formulae  vAich  occurs  at  the  internal  eigenfrequencies  of 
the  geometry  under  consideration.  Also,  a representation  of  the  most  singular 
component  in  the  Burton  and  Miller  formulation  is  developed  resulting  in  an 
integral  equation  which  is  amenable  to  numerical  solutions.  A simple  numerical 
scheme  is  introduced  which  reduces  the  large  amounts  of  computer  storage  and 
time  normally  required  for  the  solution  of  similar  problems.  This  numerical 
scheme  is  then  used  to  obtain  solutions  for  the  radiated  sound  field  generated 
by  a vibrating  piston  set  in  a sphere.  The  numerical  solutions  for  the  surface 
and  far  field  sound  patterns  are  compared  with  exact  analytical  solutions  and 
deviations  of  ten  percent  at  most  are  noted.  Since  the  symmetry  of  the  sphere 
was  not  taken  advantage  of  in  these  computations,  the  numerical  schemes  employed 
are  applicable  to  general  three  dimensional  sound  radiation  problems. 

I.  Introduction 

Hie  development  of  a simple  analytical  form  and  an  efficient  numerical  method 
for  the  prediction  of  the  characteristics  of  the  sound  fields  radiated  by  three 
dimensional  bodies  is  the  main  concern  of  this  paper.  Such  prediction  techniques 
have  a variety  of  applications  in  science  and  engineering;  for  example,  the  de- 
termination of  the  sound  fields  radiated  by  aircraft  and  underwater  vehicles. 

The  approach  developed  in  this  investigation  is  by  no  means  limited  to  acoustic 
radiation  problems  as  other  wave  phenomena  are  governed  by  similar  equations. 
Thus,  the  analytical  and  numerical  methods  employed  here  are  also  directly  ap- 
plicable to  other  fields  of  engineering  such  as  electromagnetic  antenna  theory 
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and  wave  scattering  problens. 

This  research  was  undertaken  with  the  objective  of  determining  the  applica- 
bility of  certain  integral  equation  formulations  for  the  exterior  Helmholtz 
problem  in  the  prediction  of  the  radiated  sound  fields  produced  by  three  di- 
mensional bodies.  In  principle,  integral  formulations  appear  very  attractive 
as  they  (1)  eliminate  the  need  to  consider  the  infinite  domains  normally  associ- 
ated with  radiation  problens;  (2)  reduce  the  dimensionality  of  the  problem  by 
one  (e.g,,  from  a three  dimensional  partial  differential  equation  to  a two  di- 
mensional surface  integral  equation);  and  (3)  can  readily  handle  arbitrary 
geometries  and  boundary  conditions.  All  three  properties  are  very  advantageous 
from  a computational  point  of  view  as  the  first  two  significantly  reduce  the 
computer  storage  requirement  for  solution  and  the  third  eliminates  the  need  to 
extensively  modify  the  computer  code  when  the  geometry  or  the  boundary  condi- 
tions are  altered. 

Difficulties  arise,  however,  in  the  use  of  the  Helmholtz  formulae  as  their 
solution  depends  upon  the  numerical  evaluation  of  singular,  oscillatory  inte- 
grands.^ ^ Also,  most  external  boundary  integral  representations  suffer  from  a 

nonuniqueness  of  the  solution  at  frequencies  corresponding  to  the  eigenf requen- 
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cies  of  the  associated  internal  problem  of  the  same  geometry.  Keeping  these 
difficulties  in  mind,  the  work  presented  in  this  paper  is  specifically  concerned 
with  the  following  problems:  (1)  the  development  of  an  accurate  and  efficient 
numerical  scheme  for  handling  the  oscillatory,  singular  integrands  encountered 
in  the  application  of  the  Helmholtz  formulae;  (2)  the  determination  of  the  most 
effective  procedure  for  handling  the  nonuniqueness  of  the  radiation  solution  at 
eigenvalues  of  the  associated  internal  acoustic  problem;  and  (3)  the  determination 
of  the  accuracy  of  the  resulting  solutions. 


While  there  are  many  papers  in  the  literature  (e.g.,  see  Refs.  1-8)  dealing 
with  integral  solutions  of  radiation  problems,  none  of  these  addresses  the  im- 
portant problem  of  determing  the  applicability  and  relative  efficiency  of  the 
various  integral  formulations  and  numerical  procedures  which  can  be  employed  to 
obtain  the  desired  solutions.  Instead,  most  of  these  investigations  are  limited 

} to  discussions  of  the  potential  advantages  of  the  use  of  certain  integral  form- 

\ 

ulations,  various  possible  approaches  for  the  numerical  solution  of  the  result- 
ing integral  equations,  the  nonuniqueness  of  the  solutions  of  the  integral 
equations  which  govern  the  radiation  problems,  and  potential  means  for  allevi- 
ating this  nonuniqueness  problem.  The  few  papers  (e.g.,  see  Refs.  1-6  and  9) 
that  deal  with  the  numerical  solutions  of  specific  problems  are  either  limited 
to  two  dimensional  problems  or  three  dimensional  problems  with  simple  boundary 
conditions,  such  as  perfectly  reflecting  surfaces,  which  greatly  simplify  the 
analytical  form  and  the  numerical  solution  procedure.  In  the  present  investi- 
gation the  analytical  and  numerical  schemes  are  applied  to  problems  involving 
general  boundary  conditions. 

II.  Theoretical  Considerations. 

In  this  sectioiY  an  outline  of  the  development  of  the  theory  upon  which  the 
calculations  are  based  is  presented.  The  basic  integral  representation  of  the 
solution  of  the  Helmholtz  equation  is  rigorously  developed  in  Ref.  (10)  and  will 
not  be  repeated  here;  however,  derivations  that  are  directly  related  to  the  pres- 
ent investigation  are  presented  in  detail. 

A.  Formulation  of  the  Integral  Equation 

Ihe  standard  three  dim.ensional  Helmholtz  formula  for  an  external  (radiation) 

Ki  ■ 7.10 
problem  is 

Sq  9 ^ 
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(See  Fig.l.)  where  G(P,Q)  is  a fundamental  solution  of  the  Helmholtz  equation; 
that  is: 


G(P,Q)  = 


.ikr(P,Q) 

r(P,Q) 


(2) 


and  k is  the  wave  number.  In  Eq,  (1)  represents  an  outward  normal  deriva- 
tive  with  respect  to  the  body  (i.e,  inward  with  respect  to  the  exterior  region) 
of  the  function  with  respect  to  the  variable  Q;  i.e. 


P(Q)  . Hq 

where  p is  the  acoustic  potential. 

Introducing  the  modified  admittance,  Y,  defined  as 


(3) 


Y(Q)  = MSI  / 5p(Q) 


Eq.  (1)  can  be  rewritten  as 


J J (5) 


Thus,  using  Eq.  (1)  or  Eq.  (5),  the  acoustic  potential  p(P)  at  anj  point  out- 
side the  surface  of  the  body  S can  be  determined  if  the  acoustic  potential  on 
the  surface  of  the  body  p(Q)  and  either  its  normal  derivative  (the  acous- 

tic  velocity)  or  the  admittance  Y(Q)  on  the  surface  of  the  body  are  known. 

If  the  point  P is  allowed  to  approach  the  surface  of  the  body,  Eq,  (5)  be- 
comes 


cp(Q) 


aG(P,Q) 


G(P,Q)  I ds 

an„  J q 


2 TT  o(P) 


(6) 
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if  the  surface  S is  sufficiently  smooth.  Using  Eq.  (4),  Eq.  (6)  becomes 

''  J } ‘*®q  ' (?) 

The  integral  Eqs.  (6)  or  (7)  can  now  be  solved  for  the  acoustic  potential  on 
the  surface  if  either  the  acoustic  velocity  or  the  admittance  is  known  at  each 
point  on  the  body.  Also  if  the  acoustic  velocity  is  known  over  part  of  the 
body  (i.e.  the  driving  surface)  and  the  admittance  over  the  remainder,  Eq.  (6) 
may  be  applied  on  the  driving  surface  and  Eq.  (7)  over  the  rest  of  the  body. 

Both  G and  its  first  normal  derivative  with  respect  to  the  variable  Q, 
which  appear  in  the  kernels  of  Eqs.  (6)  and  (7),  become  singular  when  the 
point  Q approaches  the  point  P on  the  surface  (See  Eq,  (2)J.  It  can  be  shown, 
however,  that  the  integrals  are  regular  in  spite  of  this  singularity  of  the 
kernels,  and  no  analytical  problems  arise  because  of  it.  However,  the  singu- 
lar kernels  do  present  numerical  difficulties  which  will  be  discussed  in 
Section  III. 

An  analytical  problem  does  arise  in  the  solution  of  Eqs.  (6)  and  (7)  when 

the  wave  number  k,  which  appears  in  the  simple  source  solution  G (See  Eq.(2).) 

approaches  a resonant  frequency  (i.e.  an  eigenvalue)  of  the  related  internal 
7 11 

problem,  ’ At  these  frequencies  Eqs.  (6)  and  (7)  do  not  yield  a unique  solution. 

B.  The  Uniqueness  Problem 

Since  the  uniqueness  problem  occurs  only  at  certain  wave  numbers  corres- 
ponding to  internal  eigenvalues  it  night  be  suggested  that  the  •'roblen  be 
simply  avoided  by  considering  only  wave  numbers  which  are  not  close  to  internal 

* It  will  be  noted  here  that  Eq.  (7)  yields  a homogeneous  set  of  equations  if 
only  the  admittance  is  kno^m.  Thus  to  obtain  a unique  solution  the  acoustic 
potential  must  be  known  on  part  of  the  body. 
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eigenvalues.  This  is  not  feasible,  however,  because:  (1)  if  the  body  is  truly 
arbitrary  in  shape  the  internal  eigenvalues  are  not  known  a priori  and  the 
corresponding  internal  problem  would  also  have  to  be  solved  in  order  to  de- 
termine what  wave  numbers  to  avoid;  (2)  the  integral  equation  is  discretized 
for  numerical  integration,  which  results  in  a system  of  algebraic  equations, 
so  that  there  is  no  longer  a specific  value  but  a range  of  values  at  which 
the  coefficient  matrix  becomes  ill-conditioned  which  results  in  large  numeri- 
cal  errors;  and  (3)  the  interval  between  successive  eigenvalues  decreases 
with  increasing  wave  number  and  it  becomes  impossible  to  stay  "sufficiently" 
far  away  from  the  internal  eigenvalues  at  high  wave  numbers  (e.g.,  k on  the 
order  of  10), 

It  has  been  suggested^  that  one  method  to  assure  the  uniqueness  of  the 
solution  is  to  obtain  an  overdetermined  system  of  algebraic  equations  by  com- 
bining the  system  of  algebraic  equations  generated  from  the  standard  integral 
equation  (e,g,,  Eq,  (6))  with  additional  algebraic  equations  generated  from 
the  integral  relation 


f 

J i if«) 

S 

q 


5GCP.Q) 


G(P,Q)  I dS  = 0 

3nq  ; q 


(8) 


where  the  point  P lies  inside  the  surface  S.  There  are  two  problems  with  this 
approach.  The  first  is  determining  the  number  of  extra  relations  required  to 
"pick-out"  the  proper  solution  from  the  set  of  possible  solutions  of  the  non- 
unique integral  equation;  and,  the  second  is  choosing  the  placement  of  the 


* If  the  admittance  is  non-zero  the  internal  eigenvalues  of  the  problem  are  in 
general  complex.  However,  even  if  only  real  wave  numbers  are  considered  the 
nonuniqueness  problem  still  exists  if  the  imaginary  part  of  the  complex  eigen- 
value is  sufficiently  close  to  zero. 
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points  which  are  used  to  generate  the  extra  relations.  As  there  is  no  known 

procedure  for  choosing  either  the  optimum  number  of  extra  relations  or  the 

points  from  which  they  are  generated,  this  method  can  not  be  relied  upon  to 

give  consistently  good  results. 
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Ursell  has  suggested  that  the  uniqueness  problem  be  avoided  by  the  use 
of  a different  fundamental  solution  (i.e.  a different  G function;  see  Eq.  (2)). 


Although  the  use  of  a different  fundamental  solution  does  not  change  the  re- 
sulting integral  equations  and  analytically  eliminates  the  uniqueness  problem 
rather  elegantly,  the  function  itself  is  difficult  to  construct  numerically  as 
it  entails  the  computation  of  infinite  series.  Thus  the  elegance  of  the  method 
is  offset  by  the  large  increases  in  computer  time  and  storage  required  for  its 
implementation,  especially  when  considering  three  dimensional  problems. 

Another  method  for  overcoming  the  uniqueness  problem  is  based  upon  the 
fact  that  a unique  solution  can  be  obtained  by  solving  a modified  integral 
equation  consisting  of  the  original  integral  equation  (6)  and  its  differen- 
tiated form^,  that  is 

I 1 1 Sjran  .2n  ^ 

J J L an  an  an  an  J q an 


p q 


(9) 


Equation  (9)  also  describes  the  behavior  of  the  acoustic  potential  on  the  sur- 
face of  the  body,  and  it  has  a set  of  related  internal  eigenvalues  which  is 
mutually  exclusive  of  the  set  of  related  internal  eigenvalues  of  Eq.  (6).  Thus, 
neither  equation  ever  fails  to  yield  a unique  solution  at  the  same  k value  as 
the  other.  Using  this  fact,  the  following  linear  combination  of  Eqs.  (6)  and  (9) 


o 


s 

q 

+ „ r r f ( |!g^  . 2£(i2X  MSI  } ,s 

J J L ^ 9n  3n  Sn  5n  J c 


P q 


= 2tt  (cp(P)  + a ) 

on 

P 


(10) 


where  a is  a coupling  constant,  should  yield  a unique  solution  for  all  values 
of  the  wave  number  k. 
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Specifically,  Burton  and  Miller  , have  shown  that  the  following  relation- 
ships exist  between  the  coupling  constant  a and  the  wave  number  k 

Im  (o)  0 -•  k real  or  imaginary 

Im  (a)  = 0 — k complex  (11) 


which  assures  a unique  solution.  Unfortunately,  the  differentiated  form  of  the 
integral  equation  (9)  contains  the  following  term 


j j 

S 

q 


fgf  dS 

an  on  q 

p q 


(12) 


which  is  strongly  singular  as  the  point  Q approaches  the  point  P.  Because  of 

its  singular  form  this  term  cannot  be  directly  integrated  numerically. 

Two  methods  for  approaching  this  problem  have  been  suggested.  The  first 

solution  is  to  "regularize"  the  singular  component  by  an  integration  of  the 

14 

entire  equation  (See  Eq.  (10).).  This  method  requires  an  excessive  amount 
of  computing  time  as  an  additional  integration  must  be  performed  over  the 
surface  of  the  body.  The  other  approach  suggests  the  use  of  a transformation 
to  interpret  the  singular  integral^^  (See  Eq,  (12).).  In  Ref.  (15)  two  alternate 
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forms  of  the  singular  integral  are  put  forth.  The  first  requires  further  manipu- 
lation to  be  of  use  as  it  contains  yet  another  singular  integral.  The  second 
requires  an  excessive  amount  of  computer  storage  space  as  it  necessitates  addi- 
tional information  that  will  allow  the  computation  of  the  tangential  derivative 
of  the  acoustic  potential  on  the  two  dimensional  surface  of  the  body.  It  must 
also  be  noted  that  the  acoustic  potential  is  the  unknown  in  most  problems  so 
that  some  differencing  procedure  is  required  to  generate  the  solution  p(Q). 


C.  Treatment  of  the  Singular  Integral 
In  this  section  the  first  relationship  developed  in  Ref.  (15)  (See  pp, 
1283-1284.)  is  used  as  a starting  point  for  deriving  the  desired  expressions. 
It  is  shown  in  Ref.  (15)  that 

I !'  1%^  ,s 

J J 5n  Sn  q 

c P q ^ 


= I ! cp(Q)  (n  •«  ) V . V G(P,Q)  dS 

J j p q p q q 

S 

q 

+ J I P(Q)  (Hp  X n^)  . (Vp  X 7^  G(P,Q))  dS^ 
S 

q 


(13) 


- i I c?(Q)  n . V X (n  x v G(P,Q))  dS 

,)  j q q p p q 

S 

q 


The  first  two  integrals  on  the  right  hand  side  are  regular;  however,  the  third 
is  not.  It  is  also  shown  that  after  some  manipulation  an  alternate  form  of  the 
third  term  is 


J J K ^q  ] • L"p  ^p  J 

S 

q 


(14) 
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This  integral  is  regular  so  that  the  singular  integral  has  been  shown  to 
be  equivalent  to  the  sum  of  three  regular  integrals.  It  should  be  noted  that 
the  first  term  in  this  integral,  |^n^  x cp(Q)J,  is  the  tangential  derivative 
of  the  acoustic  potential  on  the  surface  of  the  body  alluded  to  in  the  previous 
subsection. 

An  interesting  property  of  this  integral,  Eq.  (14),  is  that  if  the  acoustic 
potential  cp(Q)  is  a constant  on  the  surface  of  the  body,  the  integral  is  zero 
as  in  this  case  n^  x cp(Q)  = 0.  Since  the  two  formulations  are  equivalent  it 
follows  that 


- J J "q  • ^q  ^ ^"p  ^ ^p 

S 

q 

= I I ^ • ["p  ^ j 

s 

q 


(15) 


Writing  the  third  integral  on  the  right  hand  side  of  Eq,  (13)  as 


- r 1 (cp(Q)  - cp(P))  n , X (n  x v G(P,Q))  dS 
J J q q p p q 

s 

q 

(16) 

- J J "q  • ^q  ^ ^p 

S 

q 


where  c?(P)  is  a constant  with  respect  to  the  variable  Q,  we  see  that  the  last 
integral  is  identically  zero  by  setting  <p(Q)  s 1 in  Eq.  (15). 

H2nce,  the  first  term  in  Eq,  (16)  is  not  only  regular  but  it  can  also  be 
readily  integrated  numerically.  As  point  Q approaches  point  P the  entire  in- 
tegral goes  to  zero.  Thus  it  has  been  shown  that  the  singular  integral  which 
appears  in  the  "unique"  formulation  of  this  problem  can  be  expressed  in  the 
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folloving  form 


J 1 

S P ‘I 

q 


dS 


= J j cp(Q)  (Up  . n^)  Vp  . Vq  G(P,Q) 


dS 


(17) 


+ J J 9(Q)  (Hp  X n^)  . (Vp  X G(P,Q))  dS 


- J J [cp(Q)  - cp(P)]  Hq  . Vq  X (Hp  X Vp  G(P.Q))  dS^ 
S 

q 


which  will  be  used  in  the  numerical  computations  of  this  paper. 


D.  Computational  Considerations 

Because  of  the  special  form  of  the  fundamental  solution  of  the  Helmholtz 
equation,  G(P,Q)  = G(Q,P),  (See  Eq.  (2).)  certain  simplifications  can  be  made; 
specifically 


. Vq  G(P,Q)  = - Vq  G(P,Q)  = G(P,Q) 

Vp  X V^  G(P,Q)  = . Vq  X Vq  G(P,Q)  = 0 


(18) 


Using  the  above  relationships  Eq,  (17)  can  be  rewritten  as  follows: 


cn  on 

p q 


dS 

q 


= - I j cp(Q)  (np  , n^)(ik)^  G(P,Q)  dS^ 
S 

q 


" J I "q  • ^q  % G(P,Q))  dS^ 

s 

q 


(19) 


To  reduce  Eq.  (19)  to  a form  nore  amenable  to  numerical  computation  it  is 
convenient  to  let  c;(Q)  = 1 so  that 


IJ 


^0^  ds 

on  Sn  q 

p q 


J J % • G(P,Q)  dS^ 


(20) 


Using  Eq,  (20)  the  left  hand  side  of  Eq,  (19)  can  be  rewritten 


as 


jj»«) 


p q 


■ J J - "‘’'O  <“= 


p q 


(21) 


- cp(P) 


(Hp.  n^  )(ik)^  G(P,Q)  dS^ 


If  Eq.  (2)  is  employed  and  the  indicated  differentiations  are  performed 
the  right  hand  side,  Eq.  (21)  can  be  rewritten  in  the  following  form 

J J bw*  - »<^>]  - 7^,  * ^ ] 


aiCP,.q)  5r(P,Q)  _ • "q  p 

9n  an 

P q 


[r(P,Q)] 

rCP^Q)'^  [ Ik  - ] } dS^ 


, , ^ikr(P,Q) 

J J "r(P,Q)  ("p  • -q)  dSq 


(22) 


„h„.  _ - 


Using  the  results  developed  in  this  section 
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the  formulation  of  Burton  and  Miller  (See  Eq.  (10).)  reduces  to 


ikr(P,Q) 


r r (ik  . _i ^ arCP,.q} 

J J r(P,Q)  r(P,Q)  ^ Qn 

c ^ 


- " J J 


ikr(P,Q) 

(ik)  (n  . n ) dS 


r(P,Q) 


p q q 


. « j J [,(Q)  - »(p)  ] { [(iw^  - ^ 

q 


5^^(P.Q)  Sr(P..Q)  _ '"p  • \ 

dn^  an^  r(P,Q) 


(ik 


r(P,Q)  ^ } 


dS 


-IJ 


r dS 

an^  r(P,Q)  “'’q 


c^JJ  4^  (Ik-rT^)  dS 


an„  r(P,Q) 


r(P,Q)  " drip  q 


2 TT(cp(P)  + a 

P 


(23) 


The  above  equation,  although  it  may  appear  more  complicated,  is  actually 
considerably  simpler  from  a numerical  point  of  view  than  solving  Eq.  (10)  withEqs. 
ilS)  and  (14).  In  summary,  the  above  formulation  of  the  sound  radiation  problem 
provides  unique  solutions  at  all  wave  numbers  k and  contains  no  singular  integrals. 

III.  Numerical  Considerations, 

To  determine  the  radiated  sound  field  generated  by  an  arbitrarily  shaped 
three  dimensional  body,  Eq.  (23)  must  first  be  solved  for  the  distribution  of 
the  acoustic  potential  on  the  surface  of  the  body,  q)(Q).  Then  this  data  needs 
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to  be  substituted  into  Eq.  (1)  to  determine  the  radiated  sound  field.  Inspection 

of  Eq.  (23)  indicates  that  all  of  the  integrands  appearing  in  this  equation  are 

ikr(P  Q) 

both  oscillatory  and  singular  due  to  the  factor  e ‘ ' which  appears  in 

r(P,Q) 

each.  Therefore  care  must  be  exercised  in  the  numerical  representation  of  these 
kernels. 

When  considering  the  numerical  evaluation  of  an  integral  on  an  arbitrary 

two  dimensional  surface,  such  elegant  computational  methods  as  Gaussian  quad- 

rature^^’^^  (which  has  been  found  by  the  authors  of  this  paper  to  yield  accurate 
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results  in  two  dimensional  sound  radiation  problems  ) cannot  be  used  in  the 
numerical  representation  of  the  kernels.  The  only  simple  approach  available  to 
obtain  a more  accurate  representation  of  the  kernels  is  to  evaluate  them  at 
more  points  on  the  surface  of  the  body.  Unfortunately,  this  is  usually  accom- 
panied by  an  attendant  increase  in  the  size  of  the  coefficient  matrix  which 
must  be  solved  to  obtain  the  acoustic  potential.  Hie  computer  time  required 
to  solve  this  matrix  goes  up  roughly  as  the  square  of  the  number  of  unknowns 
for  most  methods  of  solution  (e,g,,  Gauss-Jordan  reduction). 

Two  considerations  enter  into  the  determination  of  the  size  of  the  coeffi- 
cient matrix;  (1)  the  heuristic  determination  of  the  number  of  points  required 
on  the  surface  of  the  body  to  represent  the  acoustic  potential  to  the  desired 
accuracy;  and  (2)  the  computer  time  and  storage  space  available  to  solve  the 
coefficient  matrix  resulting  from  the  discretization  of  the  integral  equation. 

The  storage  space  available  is  usually  much  smaller  than  the  number  of 

points  required  for  the  accurate  evaluation  of  the  singular,  oscillatory  kernels. 

in  view  of  the  above  considerations  the  following  scheme  is  used  to  obtain 
a numerical  solution.  First,  the  surface  of  the  body  is  divided  into  a number 
of  area  elements  which  corresponds  to  the  number  of  points  where  the  acoustic 
potential  is  to  be  calculated  on  the  surface.  It  has  been  determined  that 
better  results  are  obtained  in  general  if  the  area  elements  are  "regular"  (i.e. 
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not  too  elongated  in  any  direction),  although  the  exact  shape  is  unimportant, 
and  they  should  be  of  roughly  equal  area.  A point  is  then  chosen  in  the  "center" 
of  each  area  element  (usually  the  centroid  of  the  plane  figure  projected  to  the 
body  surface).  These  points  will  be  denoted  as  calculational  points  (i.e,  P 
points)  as  this  is  where  the  acoustic  potential  will  be  calculated.  Next,  each 
of  the  original  area  elements  is  subdivided  into  a number  of  smaller  area  ele- 
ments, the  sum  of  which  corresponds  to  the  number  of  points  where  the  singular, 
oscillatory  kernels  must  be  evaluated  on  the  surface  of  the  body  to  assure  their 
accurate  representation.  A point  is  thus  chosen  on  the  surface  of  the  body  in 
the  "center"  of  each  of  the  smaller  area  elements  as  before.  These  points  will 
be  denoted  as  computational  points  (i.e.  Q points)  since  the  kernel  functions 
are  computed  there. 

The  calculational  points  may  or  may  not  be  a subset  of  the  computational 
points  on  the  surface  of  the  body.  If  they  are  a subset  some  computer  space 
may  be  saved;  however,  the  computational  point  must  be  avoided  when  it  corres- 
ponds to  a calculational  point  (i.e.  when  the  point  Q corresponds  to  the  point 
P)  since  the  kernels  are  then  singular.  Thus  each  term  in  the  coefficient  matrix 
is  now  the  sum  of  a number  of  terms  generated  by  a number  of  evaluations  of  each 
kernel  function. 

Since  the  integrals  are  all  regular  a better  approximation  may  be  obtained 
by  placing  computational  points  closer  to  the  calculational  point  when  one  is 
close  to  the  singularity  of  the  kernel  function.  Thus  the  computational  area 
elements  may  be  further  subdivided  to  obtain  a more  accurate  representation  of 
the  integral  about  the  point  P. 

An  illustration  of  how  the  above  procedure  is  accomplished  is  presented 
below  using  Eq,  (23)  which  is  rewritten  in  the  following  more  compact  form 
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J j cp(Q)  A(P,Q)  dS^  + g(P)  j J B(P,Q)  dS^ 

S S 

q q 

+ J J [cp(Q)  - cp(P)]  C(P,Q)  dS^  - 2 nc?(P)  (24) 

s 

q 

- 2 + 1 J 

^ s 

q 


where  the  proper  form  of  the  integrands  can  be  readily  obtained.  Next,  Eq,  (24) 
is  discretized  as  follows: 


n n 

E cp(Q)  A(P,Q)  AS^  + p(P)  E B(P,Q)  hS 
Q=1  ^ Q=1  ^ 

QjiP 


Q= 
Q/P 


+ E lcp(Q)  - cp(P)]  C(P,Q)  hS 
Q=1  ^ -J 


(25) 


m p. 

+ cpCP)  E A(P,q)  + B(P,q)  I AS  -2^cp(P) 
q=l  ^ ^ 


= 2 D(P,Q)  A S_  + E D(P,q)  A + 2 tt  a 
Q=1  ^ q=l  ^ “^p 

Qj^P 

P = 1,2, ...,N 


where  N is  the  number  of  calculational  points;  n is  the  number  of  computational 
points  (not  including  the  subdivided  element  about  ?),  and  a is  the  number  of 
computational  points  in  the  subdivided  element.  In  the  above  representation  the 
normal  component  of  the  acoustic  velocity,  do(Q)  , on  the  body  is  assumed  known 
and  is  therefore  included  in  the  integrand  D(P,Q).  Additional  input  data  required 
to  obtain  a solution  include:  (1)  the  coordinates  of  each  computational  point; 
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(2)  the  area  associated  with  each  computational  point  and  (3)  the  outward 
normal  vector  at  each  computational  px)int.  The  information  required  at  the 
calculational  points  is  included  in  the  above. 

If  the  resulting  matrix  of  coefficients  is  large,  there  are  many  iterative 
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schemes  which  can  be  employed  in  its  solution  ; however,  if  the  matrix  is  small 
Gaussian  elimination  with  back  substitution  may  be  used.  Once  the  acoustic  po- 
tential is  determined  on  the  surface,  Eq.  (1)  may  be  used  to  generate  the 
acoustic  potential  at  any  point  in  the  field  surrounding  the  body.  In  this  com- 
putation the  point  Q never  coincides  with  point  P and  the  integrands  are  never 
singular;  however,  they  are  still  oscillatory  and  care  must  still  be  taken  to 
get  an  accurate  representation  of  the  integrals. 

Due  to  the  availability  of  analytical  solutions  for  comparison  purposes, 
the  developed  numerical  procedure  has  been  applied  to  predict  the  sound  radi- 
ated by  a sphere.  However,  it  should  be  reiterated  that  no  advantage  \-as  taken 
of  the  sphere's  relatively  simple  geometry  (i.e.  its  symmetry)  in  the  numerical 
computations.  Once  the  needed  input  data  was  generated  it  was  treated  like  any 
other  arbitrarily  shaped  three  dimensional  body.  The  sphere  was  subdivided  into 
80  triangles  by  first  taking  an  icosahedron  (a  three  dimensional  figure  whose 
surface  consists  of  20  equilateral  triangles)  inscribed  in  a unit  sphere  and 
dividing  each  triangle  into  four  others  (See  Ref.  (1),  pp.  1630-1631.).  This 
was  accomplished  by  finding  the  midpoint  of  each  side  of  each  triangle  and  pro- 
jecting it  to  the  surface  of  the  sphere  as  shown  in  Fig,  2.  The  centroid  of 
each  triangle  was  then  found  and  also  projected  to  the  surface  of  the  sphere. 
These  80  points  correspond  to  the  previously  described  calculational  points. 

To  obtain  the  computational  points  this  method  was  simply  repeated  three 
more  times  yielding  5120  points.  Around  each  calculational  point  the  three  sur- 
rounding triangles  were  then  divided  once  more,  and  at  the  calculational  point 
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itself  the  triangle  was  divided  twice  more  as  shown  in  Fig.  3.  The  spherical 
area  was  then  computed  for  each  computational  triangle  (i.e,  the  sum  of  the 
areas  of  all  the  computational  triangles  is  4 n,  the  surface  area  of  the  unit 
sphere).  This  yielded  all  the  geometrical  input  data  required,  as;  for  a unit 
sphere  the  coordinates  of  the  computational  pioints  and  the  elements  of  the 
out\7ard  normals  are  the  same  in  rectangular  coordinates. 

IV, Results  and  Discussion, 

In  the  calculations  performed  in  this  study  the  surface  of  the  radiating 
sphere  was  divided  into  two  parts.  On  one  part  (the  driving  surface)  the  normal 
acoustic  velocity,  ^ , was  specified  while  on  the  other  part  (the  admittance 
surface)  the  modified  admittance  function  Y,  defined  by  Eq.  (4),  was  specified 
indicating  either  sound  absorption  or  amplification  by  this  part  of  the  surface 
(See  Fig.  4.).  The  sphere  was  chosen  for  this  study  as  exact  analytical  solutions 
can  be  obtained  for  comparison  with  the  numerical  solutions  obtained  by  solving 
the  integral  equations. 

Using  the  well  known  separation  of  variables  technique  it  can  be  shown  that 
the  acoustic  potential  for  the  sphere  can  be  represented  as  follows 


{(iV)'"  ^ F„(l)  } {3°:  : 


(26) 


either  on  or  in  the  field  surrounding  the  surface  of  the  sphere.  In  the  above 

expression  T\  = cos  9,  ^ = kr,  h^  is  a spherical  Haokel  function  of  order  m, 

and  is  a Legendre  polynomial  of  degree  n.  It  should  be  noted  that  when 

n = 0 all  5 dependence  drops  out  so  that  the  problem  becomes  axi-symmetric , 

20 

It  can  also  be  shown  that  the  acoustic  potential  for  a piston  vibrating 

in  an  otherwise  hard  (i.e,  Y = 0)  sphere  is  given  by 


<p(r,0.§) 


l^/2k)  - 4+1  (4>  j{ 

m=U 


4* 


(27) 
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both  on  the  surface  of  the  sphere  and  in  the  field  surrounding  it.  In  Eq,  (27), 

0^  denotes  the  edge  of  the  piston  set  in  the  sphere,  a is  the  radius  of  the 
sphere,  = ka,  = cos  9^,  P ^ ^ (when  n = 0),  and  the  remaining 

quantities  are  the  same  as  those  appearing  in  Eq.  (26).  The  solution  is  al- 
ways axi-syimnetric  as  there  is  no  ^ dependence.  Also,  both  solutions  (i.e, 

Eqs,  (26)  and  (27))  represent  radiated  sound  fields  as  they  satisfy  the 
necessary  radiation  conditions  when  r -•  =>. 

In  all  the  calculations  performed  in  this  study  the  acoustic  velocity 
is  specified  on  a quarter  of  the  sphere's  surface  as  shown  in  Fig,  4,  Also, 
the  sphere  is  of  unit  radius  (i.e,  a = 1)  and  the  coupling  constant  a,  re- 
quired in  applying  the  method  of  Burton  and  Miller,  has  been  taken  as  the 
pure  imaginary  number  i as  k is  a real  number  (See  Eq.  (11).). 

The  radiative  fields  computed  in  this  study  are  sunmarized  in  Table  I 
where  the  assigned  values  of  m and  n describe  a specific  solution  (See  Eq,  (26).). 
For  each  investigated  case  Table  I contains  the  exact  solutions  on  the  surface 
and  far  fields,  the  input  boundary  conditions  derived  from  the  known  exact 
solution,  and  the  average  percent  error  obtained  by  comparing  the  exact  and 
computed  solutions. 

To  check  the  numerical  approach  and  computer  code  Case  ir  1 was  investigated 
initially  (See  Table  I.).  Under  these  conditions  the  analytical  solution  for  the 
surface  potential,  tp  = h^(kr),  is  a constant  as  is  the  input  data  (i.e,  ^ and  Y). 
A comparison  between  the  numerical  and  exact  solutions  for  the  amplitude  {cp|  of 

'is 

the  acoustic  potential  on  the  surface  of  the  sphere  and  in  the  far  field  is 
presented  in  Fig,  (5)  where  excellent  agreement  between  the  two  solutions  is 
noted.  The  computation  of  the  surface  solution  required  three  minutes  of  computing 


* The  amplitude  was  chosen  for  comparison  as  both  the  real  and  imaginary  parts 
of  the  acoustic  potential  show  similar  trends  and  errors. 


time  on  the  Georgia  Tech  CDC  Cyber  70  model  74  computer.  This  tine  is  indicative 
of  all  the  cases  considered  in  this  study.  Ihe  far  field  distribution  of  the 
acoustic  potential  was  calculated  using  both  the  exact  surface  distribution 
and  the  calculated  surface  distribution.  The  far  field  calculation  using  the 
exact  surface  distribution  was  done  as  a check  on  the  computer  code  and  the 
results  agreed  with  the  exact  solution  (obtained  from  Eq,  (26))  to  seven  sig- 
nificant figures  in  both  the  real  and  imaginary  parts.  To  calculate  both  distri- 
butions simultaneously,  under  two  minutes  of  computing  time  was  required  which 
is  also  indicative  of  all  the  cases  run. 

The  second  investigated  solution  (See  Case  # 2,  Table  I.)  was  the  same 
solution  as  Case  # 1 but  with  a wave  number  k = tt  which  coincides  with  the  first 
internal  eigenvalue  of  the  sphere.  This  case  was  run  to  check  the  validity  of 
the  theory.  Again  the  solution  for  the  surface  potential  and  the  input  data  are 
constants  both  on  and  off  the  surface  of  the  sphere.  A comparison  between  the 
computed  and  exact  solutions  on  the  sphere  surface  and  in  the  far  field  are  pre- 
sented in  Fig.  (6).  Examination  of  this  figure  indicates  that  in  this  case  the 
agreement  is  not  as  good  as  in  Case  ir  I , although  the  average  error  was  still 
under  ten  percent.  Ihe  far  field  distribution  of  the  acoustic  potential  was  cal- 
culated as  before.  Employing  the  exact  surface  distribution  the  results  compared 
with  the  exact  far  field  solution  to  four  significant  figures  in  both  the  real 
and  imaginary  parts.  In  examining  Fig.  6 it  is  interesting  to  note  that  there 
was  no  increase  in  the  error  from  the  surface  distribution  to  the  far  field.  It 
was  found  that  there  was  no  significant  increase  in  error  from  the  surface  to 
Liie  far  field  distribution  in  any  of  the  cases  run. 

The  exact  solution  for  Case  # 3 (obtained  from  Eq.  (26)  with  the  data  in 
Table  I)  is  = hj^  (kr)  cos  0.  In  contrast  to  the  previous  cases  this  solution 
is  0 dependent.  A comparision  between  the  exact  and  numerical  solutions  on  the 
surface  of  the  sphere  is  presented  in  Fig.  7,  and  a far  field  comparison  is  pre- 
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sented  in  Fig,  8.  Examination  of  the  data  shows  that  in  this  case  the  error 
was  actually  reduced  in  going  from  the  surface  to  the  far  field.  Furthermore 
in  the  far  field,  the  calculated  surface  distribution  gives  better  results 
than  does  the  exact  surface  distribution, 

Ihe  next  case  investigated  (Case  # 4)  was  run  at  the  second  internal  eigen- 
value of  the  sphere,  k = 4.49340946.  The  results  for  the  acoustic  potential  on 
the  surface  of  the  sphere  are  presented  in  Fig.  9.  As  can  be  seen  the  results 
deteriorate  somewhat  at  an  internal  eigenvalue  of  the  problem.  It  is  interest- 
ing to  note  that  the  error  increases  with  0 and  it  reaches  its  maxinum  value 

at  0 = 180°,  the  center  of  the  admittance  surface.  The  results  for  the  far 

field  are  presented  in  Fig,  10,  All  of  the  cases  considered  so  far  were  axi- 
symraetric,  that  is,  there  was  no  § dependence;  a property  that  was  also  retained 
by  the  developed  numerical  solutions. 

This  next  case  (Case  # 5)  is  truly  three  dimensional  as  there  is  a § de- 
pendence in  the  solutions.  Referring  to  Eq.  (26)  and  Table  I the  exact  solu- 
tion on  the  surface  of  the  sphere  is  found  to  be  cp  = hj^(kr)  sin  Q sin  §.  The 

average  percent  error  was  not  calculated  in  this  case  due  to  the  zeros  which 
appear  in  the  exact  solution,  but  the  errors  remained  small  (i.e.  under  ten 
percent).  The  far  field  distribution  of  the  acoustic  potential  was  calculated 
and  no  increase  in  error  was  detected. 

In  the  next  three  cases  a hard  sphere  (i.e.  Y = 0 on  the  admittance  sur- 
face) with  a unit  driver  (i.e.  ^ = 1 on  the  driving  surface)  was  considered, 

on 

The  exact  solutions  for  these  cases  can  be  obtained  from  Eq.  (27).  In  these 
studies  the  wave  number  k was  varied  to  determine  the  value  of  k at  which  the 
accuracy  of  the  solution  deteriorates  for  a fixed  number  of  80  calculational 
points  that  was  used  in  these  numerical  studies. 
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Tlie  solution  for  Case  # 6 where  k = 2,  are  presented  in  Fig,  11  for  both 
the  surface  distribution  of  |cp  j and  the  far  field  solution.  In  this  case  the 
far  field  is  considered  to  be  at  kr  = 100, 

In  this  next  case  (Case  # 7)  the  wave  number  is  increased  to  k = 5.  As 
can  be  seen  from  Fig,  12,  the  error  is  still  under  ten  percent  on  both  the 
surface  of  the  sphere  and  in  the  far  field  (i.e,  kr  = 100),  When  the  wave 
number  is  increased  to  k = 10  (Case  # 8)  the  error  becomes  rather  large.  The 
average  error  in  the  calculation  of  the  surface  potential  is  sixty  percent. 

So  it  can  be  seen  that  there  are  not  enough  calculated  points  to  accurately 
represent  the  potential  function,  cp.  The  far  field  is  calculated  at  kr  = 100 
and  the  error  drops  a bit  but  it  still  remains  hi^  at  twelve  percent. 

V.  Summary. 

In  summaiy,  a solution  approach  has  been  developed  in  this  paper  which 
may  be  used  to  yield  a unique  solution  for  the  distribution  of  the  acoustic 
potential  on  the  surface  of  an  arbitrary  three  dimensional  body  at  all  values 
of  the  wave  number.  Also,  a numerical  scheme  was  developed  to  solve  the  equa- 
tion accurately  and  efficiently.  Computer  programs  were  run  to  verify  the 
applicability  of  the  developed  solution  method  and  to  find  its  limit  of  accu- 
racy for  a fixed  number  of  points.  The  procedure  was  found  to  be  both  accurate 
and  versatile  as  the  computer  code  required  no  major  modifications  to  handle 
the  various  boundary  conditions  imposed  on  the  surface  of  the  body. 
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Figure  2.  Method  of  Dividing  the  Surface  of  the  Sphere 
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Figure  3,  Division  of  the  Surface  of  the  Sphere  Around  a 
Calculational  Point, 
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Figure  8,  Exact  and  Calculated  Values  of  | cp  | using  the 
Exact  and  Calculated  Surface  Distribution  of  p 
for  p = hi(2r)  cos  0 in  the  Far  Field  of  the 
Sphere, 
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Figure  9.  Exact  and  Calculated  Values  of  ) p | for  y = hj(4.49r)  cos  c? 
on  the  Surface  of  the  Sphere. 
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Figure  10.  Exact  and  Calculated  Values  of  | cp  | using  the 
Exact  and  Calculated  Surface  Distribution  of  cp 
for  cp  = h2^(4.49r)  cos  3 in  the  Far  Field  of  the 
Sphere. 
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Figure  11.  Exact  and  Calculated  Values  of  | ® | for  a Hard  Sphere 
with  a Unit  Driver  (k  = 2)  on  the  Surface  and  in  the 
Far  Field. 
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Figure  12,  Exact  and  Calculated  Values  of  | cp  | for  a Hard  Sphere 
with  a Unit  Driver  (k  = 5)  on  the  Surface  and  in  the 
Far  Field. 
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Abstract 

A general  analytical  method  for  determining  the 
radiated  sound  fields  fr^mi  axisymmetric  surfaces  of 
arbitrary  cross  section  with  general  boundary  con- 
ditions is  developed,  Ihc  method  is  based  on  an  in- 
tev.ral  representation  cf  the  external  solutions  of 
the  Helmholtz  equation  valid  at  all  wave  numbers, 

Ihe  axi  svTBTie  t r ic  formulation  of  the  problem  reduces 
Its  solution  to  the  numerical  evaluation  of  line 
integrals  by  Gaussian  quadrature.  The  applicability 
ot  ihi  solution  approach  for  both  a sphere  and 
finite  cylinder  is  demonstrated  by  comparing  the 
nuru  rival  results  with  exact  analytical  solutions 
for  both  discontinuous  and  ci^ntinuous  boundary  con- 
ditions, 

I . Introduction 

To  reduce  tlie  noise  radiated  to  the  community 
from  turbofan  inlets,  the  effects  of  sound  suppres- 
i.ion  material  in  the  inlet  and  the  spatial  distri- 
l>ution  of  tlie  sound  source  on  the  radiated  sound 
levi-ls  and  patterns  must  he  determined.  Analytical 
techniques  for  predicting  these  effects  must  be 
capable  of  dealing  with  general  axisymmetric  geome- 
trit's  and  k^ompltcated  boundary  conditions  which  are 
encountered  in  mu i t iply- I ined  inlets.  To  determine 
the  radiated  sound  field,  an  additional  requirement 
is  that  the  methods  be  applicable  to  infinite  do- 
mains, Elu'  objective  of  this  paper  is  to  develop  a 
gem  ral  analytical  method  for  determining  the  radi- 
ated sound  fields  from  axisymmetric  surfaces  of 
arbitrary  cross  section  and  with  general  boundary 
conditions, 

Ihe  method  used  in  this  investigation  is  based 
t>n  an  integral  form  of  the  solutions  of  the  Helm- 
holtz equation, With  this  formulation  the  acous- 
tic potential  anywliere  external  to  the  surface  can 
b(‘  found  once  the  potential  distribution  on  the 
surface  is  known.  Thus,  to  determine  the  radiated 
sound  field  the  problem  reduces  to  the  determina- 
tion of  the  distribution  of  the  acoustic  potential 
on  the  two-dimensional  surface  of  the  geometry  un- 
der consideration  instead  of  solving  the  Helmholtz 
equation  in  the  surrounding  infinite  three  dimen- 
sional domain. 

It  has  been  previously  shown^”^  chat  when  ap- 
plied to  exterior  sound  radiation  problems  the  solu- 
tion technique  fails  to  produce  unique  solutions  at 
frequencies  correspionding  to  interior  eigenvalues 
of  the  geometries  under  consideration.  Unless  spe- 
cial precautions  are  taken,  stralglic-forward  numeri- 
«..il  solutions  of  the  developed  integral  equation  at 
Ircquencics  close  to  the  eigenvalues  of  the  Internal 
problem  produce  large  errors,  A technique  proposed 
by  Hurt«)n  and  Miller^  for  avoiding  this  uniqueness 
problem  and  the  associated  numerical  errors  Is  used 
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in  this  investigation.  TTiis  technique  involves  a 
reformulation  of  the  "classical"  integral  equation 
and  the  solutions  obtained  are  valid  at  all  frequen- 
cies. 

The  resulting  integral  equation  for  tlie  surface 
acoustic  potential  is  solved  numerically  and,  for 
axisymmetric  geometries,  the  equation  reduces  to 
the  evaluation  of  a line  integral.  Ihus,  the  axi- 
symmetric case  can  be  reduced  to  an  equivalent  one- 
dimensional problem,  lliis  equation  is  discretized 
and  the  resulting  system  of  algebraic  equations  is 
solved  using  complex  Causs-Jordan  elimination.  Since 
the  coefficient  matrix  involves  the  free  space 
Green's  function,  which  becomes  singular  as  two 
points  on  the  surface  approach  one  another,  numeri- 
cal techniques  are  presented  wliich  can  deal  with 
these  singularities  and  yield  accurate  results. 
Gaussian  integration  is  used  to  increase  the  accura- 
cy of  the  solution  without  significant  penalties  in 
computer  storage  and  time  requirements,  Ihe  applic- 
ability of  the  integral  formulation  and  the  accuracy 
of  the  numerical  techniques  are  demonstrated  by 
computing  the  surface  and  far  field  distributions 
of  the  acoustic  potential  on  both  a s|Vherc  and  a 
finite  cylinder.  Ihe  numerical  results  are  compared 
with  known  exact  solutions  generated  by  the  separa- 
tion of  variables  technique.  Surfaces  with  spatially 
varying  forcing  functions  and  admittances  are  con- 
sidered, for  different  tangential  modes,  Co  evalu- 
ate the  capability  of  the  integral  approach  to 
handle  boundary  conditions  of  a general  nature. 

With  the  sphere,  agreement  between  computed  and 
exact  results  is  to  three  significant  figures.  For 
Che  cylinder  agreement  is  to  two  significant  fig- 
ures, Ihe  effect  on  the  accuracy  of  discontinuous 
boundary  conditions  involving  nonzero  admittances 
over  the  surface  and  of  the  corners  encountered  in 
the  cylindrical  configuration  are  also  presented, 

II.  Theory 

In  this  section  the  general  three  dimensional 
integral  representation  of  the  solutions  of  Che 
Helmholtz  equation  is  developed  for  application  to 
radiation  problems,  Ihis  particular  formulation 
yields  unique  solutions  at  all  frequencies  and  does 
not  have  strong  singularities  wliich  are  difficult 
to  handle  numerically.  Ihe  general  integral  equa- 
tion is  then  specialized  for  axisymmetric  geome- 
tries, A more  detailed  development  is  given  in  Ref. 
5. 

General  Theory 

Beginning  with  the  three  dimensional  Helmholtz 
equation  which  governs  the  spatial  dependence  of 
the  acoustic  field  for  sinusoidal  oscillations 

7^cp  + k^(0*0  (1) 

where  is  tlie  acoustic  potential  and  k is  the  wave 
number.  The  standard  integral  representation  of  the 
exterior  solutions  is  found  to  be^*^ 
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rVio  term  represents  an  outward  normal  deriva- 
tive with  q respect  to  the  body  as  shown  in  Fig. 
1 ; that  is 
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Hic  soundest  approach  from  an  anilytical  point 
of  vit*w  ib  given  by  Burton  and  MilUr^  who  h.ivc 
.suggested  tht  use  of  the  differencial  forms  of  Kq. 
(5)  which  governs  the  spatial  dependt  lue  of  the 
acoustic  velocity. 
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Also,  G(P,Q)  is  a fundamental  three  dimensional 
solution  of  the  Helmholtz  equation  and  is  taken 
be  tlie  free  space  Green’s  Function  for  a point 
source^  defined  as 
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point  on  the  surface  of  the  body  then  the  acoustic 
potential  may  be  calculated  anywhere  in  the  exterior 
domain. 


To  solve  for  the  surface  potential,  the  point  P 
is  moved  to  the  surface  of  the  body.  Equation  (2) 
then  bectmics 
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if  the  surface  of  the  body  Is  sufficiently  smooth. 
Introducing  a modified  admittance  function  defined 
as 
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Hiis  equation  can  al.t)  be  solved  for  C5(Q)  once  the 
normal  velocity  or  admittance  is  specified  at  Che 
surface.  However,  this  equation  has  its  own  set  of 
associated  eigenvalues  at  which  unique  solutions 
cannot  be  obtained.  Burton  and  Miller  suggest  tak- 
ing a linear  combination  of  the  two  equations  to 
obtain 
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Since  the  two  sets  of  associated  internal  eigen- 
values are  mutually  exclusive  the  linear  combination 
of  equations  should  yield  unique  solutions  if  the 
complex  coupling  constant  a is  properly  chosen.  It 
is  shown  Chat  a must  meet  the  following  restrictions 
to  guarantee  chat  Kq,  (8)  yield  unique  solutions 
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Eq,  (5)  can  be  written  as 
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where 


J " J. 


If  either  the  acoustic  velocity  or  the  admittance  Is 
known  at  each  point  on  the  surface  of  the  body  then 
the  acoustic  potential  may  be  calculated  at  each 
point  using  Eq,  (7), 

Unfortunately  this  equation  does  not  yield 
unique  solutions  when  the  wave  number  k is  an  in- 
ternal eigenvalue  of  the  body  under  consideration. 
Since  the.se  eigenvalues  are  not  known  a priori  for 
general  bodies  Che  formulation  cannot  be  relied 
upon  to  give  consistently  good  results.  There  arc 
many  papers  In  the  literature^»^»^  dealing  with 
this  problem.  The  relative  merits  and  shortcomings 
of  the  methods  employed  are  discussed  in  detail  in 


A problem  arises  in  the  numerical  solution  of  Eq, 

(8)  as  the  third  term  on  the  right  hand  side  is 
strongly  singular  in  its  present  form  as  the  point 
Q approaches  the  point  P on  the  surface  of  the  body, 
Meyer,  Bell  and  Zinn^  have  shown  that  this  diffi- 
culty can  be  overcome  by  the  proper  interpretation 
of  this  singular  term.  Employing  a vector  transfor- 
mation^  and  taking  the  Cauchy  Principle  Value  Eq, 
(8)  is  shown  to  be  equivalent  to 
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All  I'i  the  Lcrms  in  Kq,  (10)  are  now  regular  and 
llurciort  arc  directly  intcgrablc;  however,  all  the 
integrands  are  oscillatory  and  singular  so  that  care 
must  be  taken  in  their  nvui\erical  approximation. 

Axisymmetric  Formulation 

Wlien  dealing  with  a body  of  revolvition  as  shown 
in  Fig,  2 an  axisyrnnetr  ic  formulation  of  the  prob- 
lem is  advantageous Ihis  being  the  case  an  element 
of  area  he  Cannes 
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where  S is  the  distance  along  the  perimeter  of  the 
surface  in  the  t-z  plane. 

Assuming  an  acoustic  vt  locity  distribution  of 
the  form 

^ = v(S)  cos  m « (12) 
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Substituting  Kqs,  (15)-(17)  into  Fq.  (14)  gives 
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Kq,  (10)  becomes 
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where  1.  is  the  length  of  the  generating  line  of  the 
surface  of  revolution.  The  S-0  coordinate  directions 
have  now  been  effectively  uncoupled  so  that  the 
problem  has  been  reduced  to  the  evaluation  of  line 
integrals  in  the  coordinate  directions  on  the  sur- 
face of  the  body.  This  formulation  does  not  restrict 
the  form  or  type  of  boundary  conditions  on  the  body; 
it  merely  assumes  that  the  boundary  conditions  can 
be  represented  by  a sum  (expanded  in  a set)  of  tan- 
gential inodes. 


III.  Re suits 
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Now,  three  sets  of  functions  arc  defined: 
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Numerical  results  have  been  obtained  for  a 
sphere  and  cylinder  using  the  numerical  technique 
described  in  Ref,  9,  Basically,  this  method  con- 
sists of  first  specifying  the  i-z  coordinates  and 
the  normal  vector  at  each  point  on  the  surface. 

From  these  quantUf^s  the  distances  r and  the  nor- 
mal derivatives  2 can  be  obtained.  The  integral 

in  Eq,  (18)  is  ^”q  then  separated  into  n inte- 
grals taken  over  subintervals  of  length  -t/n.  The 
acoustic  potential  is  assumed  constant  over  each 
subinterval  and  the  integrations  are  performed  nu- 
merically using  a four-point  Gauss-Legendre  quadra- 
ture in  the  t-z  plane,  A twenty-point  Gauss-  Legen- 
dre quadrature  formula  is  used  In  Che  circumferen- 
tial direction. 

Exact  results  were  obtained  using  separation  of 
variables,^  To  eliminate  the  need  for  evaluating 
the  resulting  infinite  series,  Che  normal  velocity 
and  admittance  distributions  were  selected  so  that 
only  one  term  in  the  series  remains. 
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twenty  points  on  the  s{)luTe,  Tlu* 
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where  r is  the  distance  from  the  origin  to  a ix>int 
P on  the  surface.  As  shown  in  Fig,  1,  with  cv  = 0 
the  computed  magnitudes  of  tlie  acoustic  potential 
are  in  error  by  1 2 per  cent  at  wave  numbers  close 
to  the  internal  eigenfrequencies  of  tt,  2r  and  3-*“ . 
lliese  results  are  those  that  would  be  obtained  from 
tq.  (3),  The  relatively  large  errors  are  expected 
from  the  analysis  of  Burton*  and  from  previous  in- 
vestigations using  Kq.  (3).^*^  Burton  proves  that 
setting  the  imaginary  part  of  cy  nonzero  guarantees 
unique  solutions  when  Kq.  (18)  is  used.  Although 
the  maximum  error  is  reduced  for  = i to  less  than 
4 per  cent  wlicn  the  nond  imens iona  I frequency  k is 
less  than  seven,  significant  errors  arc  still  evi- 
dent at  the  higher  frequencies  as  shown  in  Fig.  3, 

in  this  study  consistently  good  results  are  ob- 
tained only  when  a = i/k.  In  Fig,  i,  the  computed 
and  exact  results  with  cy  = i/k  agree  to  three  sig- 
nificant figures.  The  reason  for  this  behavior  is 
currently  under  investigation;  however,  for  all  the 
cases  presented  hereafter  this  value  of  cy  is  chosen 
and  the  exact  surface  distribution  is  given  by  Kq, 
tly)  when  m = 0. 

A problem  of  more  practical  importance  is  Che 
finite  axisymmetric  duct  since  this  surface  approxi- 
mates an  engine  configuration,  Tl^e  surface  potential 
distributions  are  presented  in  Fig.  4 for  a zero 
admittance  everywhere  on  the  surface.  The  velocity 
distribution  is  specified  over  Che  entire  surfact* 
and  the  potential  given  by  Kq,  (19)  has  a magnitude 
independent  of  frequency  and  a |>hase  linearly  pro- 
{Kjrtional  to  the  frequency.  In  Fig,  4 the  magnitude 
and  phase  arc  plotted  against  the  distance  along 
the  perimeter  S,  The  largest  errors  in  the  magni- 
tude of  Che  potential  of  about  10  percent  occur  on 
the  ends  of  the  cylinder  and  at  the  corners.  The 
results  at  the  ends  can  be  improved  without  increas- 
ing the  number  of  points  by  area  weighting  rather 
than  by  taking  equidistant  {)oints  along  the  perime- 
ter, The  errors  at  the  corners  are  caused  by  the 
discontinuous  normal  derivative  in  going  from  the 
cylinder  to  the  end.  The  errors  in  the  phase  are 
less  than  four  per  cent  in  all  cases.  The  errors  in 
magnitude  increase  with  increasing  frequency,  but 
even  when  k » 10  the  numerical  results  are  within 
10  per  cent  of  the  exact  solutions, 

in  most  practical  problems  the  boundary  condi- 
tions are  discontinuous  with  the  acoustic  velocity 
or  potential  specified  over  part  of  the  surface  and 
the  admittance  over  the  rest.  To  determine  the  ef- 
fect of  the  discontinuities  on  the  numerical  results, 
a cylinder  with  the  velocity  specified  on  the  ends 
and  tile  admittance  specified  in  the  center  was  in- 
vestigated and  the  results  are  presented  In  Fig.  5, 
Although  the  errors  of  the  numerical  results  for 
this  case  are  increased  compared  with  the  errors 
shown  in  Fig,  4,  the  errors  are  within  10  per  cent 
for  values  of  k less  than  5.  However,  when  k * lO 
errors  of  up  to  40  per  cent  in  the  magnitude  of  the 
potential  are  encountered  close  to  the  discontinu- 
ity in  the  boundary  condition.  TTiis  result  suggests 


that  more  ix>ints  need  to  be  taken  at  higher  frequen- 
cies with  discontinuous  boundary  conditions  present. 

At  higher  tangential  modes,  the  variation  in 
the  circumferential  direction  behaves  as  cos  m 
where  m = 0,1,2 To  check  the  numerical  integra- 

tion scheme  in  the  circumferential  direction,  the 
surface  acoustic  {XJtential  was  computed  for  m = 1 
and  m = 2.  Hie  results  are  presented  in  Fig.  6 for 
k = * with  the  velocity  specified  and  the  admi  t tance 
zero  everywhere  on  the  surface.  Hie  computed  and 
exact  results  are  in  agreement  to  within  two  per 
cent  for  both  m = 1 and  m = 2, 

it  has  been  shown^  that  once  tlic  surface  poten- 
tial has  bei  n accurately  computed,  tlie  far  field 
can  be  determined  to  at  least  the  accuracy  of  the 
surface  potential.  iHis  result  is  confirmed  by  the 
data  presented  in  Fig,  7 for  a cylinder  with  the 
velocity  specified  ever) where  on  the  surface  at 
k = 2,  Tlie  results  at  20  radii  from  the  surface  are 
in  agreement  with  exact  rc&ulCs  to  within  one  per 
cent  even  though  the  surfac"  ■’'rors  at  some  points 
is  above  two  per  cent.  Data  . ig,  8 show  that 
accurate  results  are  obtaineu  at  distances  greater 
than  one  integration  stepsize  from  the  surface.  At 
closer  distances  errors  from  the  numerical  evalua- 
tion of  the  singularity  in  the  Creen's  function  de- 
fined by  Kq,  (4)  leads  to  large  errors, 

I V , Summary  and  Conclusions 

An  integral  solution  of  the  Helmholtz  equation 
is  developed  for  use  in  acoustic  radiation  problems. 
Unlike  previous  formulations  which  give  poorr'^sults 
at  frequencies  corresponding  to  eigenf requencies  of 
the  surface  under  consideration,  the  formulation 
used  in  this  study  is  valid  at  all  frequencies,  'Hie 
surface  potentials  computed  numerically  for  a spliere 
and  cylinder  using  20  points  along  the  perimeter 
are  accurate  to  within  ten  per  cent  for  no.idlmen- 
sionaL  frequencies  ka  of  from  one  to  ten  where  k is 
the  wave  number  ano  a is  the  radius  of  the  sphere 
or  cylinder.  For  discontinuous  boundary  conditions, 
the  numerical  and  exact  values  are  in  agreement  to 
within  10  per  cent  for  ka  < 5.  At  higher  frequencies 
the  results  are  as  much  as  40  per  cent  in  error  at 
the  point  of  discontinuity  which  suggests  taking 
more  points  in  evaluating  the  integral  Helmholtz 
equation  to  incre  ise  the  accuracy  when  discontinu- 
ous boundary  conditions  are  specified.  At  distances 
greater  than  the  numerical  integration  stepsize, 
the  far  field  results  are  at  least  as  accurate  as 
the  corresponding  surface  potential  solutions. 
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Figure  6,  (.'omputed  Surface  Potential  for  a Cylinder 
at  the  First  and  Second  Tangential  Modes 
for  k = 2 and  20  Points, 


V 


Figure  3,  ol  Oiscont  inuous  Boundary  Conditions 

on  the  Computed  Surface  Potential  for  a Cylinder- 
2U  I’oints. 


Figure  7,  Computed  Far  Field  Potential  Distribution 
for  a Cylinder  at  k = 2,  m = 0,  and  20 
Radii  from  the  Center  of  tlie  Cylinder. 


Figure  8.  Dependence  of  the  Far  Field  Solu- 
tion of  a Cylinder  upon  the  Dis- 
tance from  the  Surface  for  k * 2 
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Predicting  the  Acoustics  of  Arbitrarily  Shaped  Bodies 
Using  an  Integral  Approach 


William  A.  Bell,*  William  L.  Meyer, t and  Ben  T.  ZinnJ 
Georgia  Institute  of  Technology,  Atlanta,  Ga. 


An  inlef>ral  Milutinn  of  the  Helmollz  equation  U developed  for  predicting  the  acoustic  properties  of  arbitraril.v 
shaped  bodies.  With  the  integral  formulation,  the  acou.stic  potentials  at  the  surface  are  solved  independently  of 
the  internal  acoustic  Held  which,  effectively,  reduces  the  dimensionality  of  the  problem  by  one.  Considerable 
reductions  in  computation  time  and  storage  requirements  are  thus  achieved.  Efficient  numerical  techniques  for 
solving  the  resulting  algebraic  equations  are  presented.  Numerical  results  obtained  for  the  two-dimensional 
problems  of  a circle  and  a rectangle  agree  In  within  one  percent  with  available  exact  solutions.  The  modes  of  a 
star-shaped  configuration  and  a duct  with  a right-angle  bend  are  also  determined  to  demonstrate  the  ap- 
plicability of  this  method  to  complicated  geometries  and  general  boundary  conditions.  Ihe  acoustic  properties 
of  a sphere  are  investigated  using  an  axisymmelric  formulation.  With  Ihe  axisymmelric  formulation  Ihe 
numerical  and  exact  results  agree  to  three  significant  figures. 


I.  Introduction 

The  prediction  of  the  acoustics  of  arbitrarily  shaped 
bodies  has  a variety  of  applications  in  aerospace 
engineering.  Among  them  are  the  determination  of  the  in- 
ternal and  radiated  sound  fields  from  airbreathing  propulsion 
systems  and  the  investigation  of  the  stability  limits  of  rocket 
combustors.  These  studies  are  concerned  with  obtaining 
solutions  to  the  Helmholtz  equation,  which  is  derived  from 
the  wave  equation  when  a sinusoidal  time  dependence  is 
assumed  and  which  describes  the  spatial  dependence  of  the 
oscillations.  This  equation  is  included  in  most  standard  texts 
on  differential  equations  of  mathematical  physics  (Ref.  1 , Ch. 
11)  and  has  been  extensively  studied  in  both  differential  and 
integral  form.  The  differential  form  is  currently  the  most 
widely  used. 

In  differential  form,  solutions  of  the  Helmholtz  equation 
can  be  obtained  by  separation  of  variables.'  * This  method 
involves  series  expansions  of  the  solutions  m terms  of 
eigenfunctions  of  the  system.  Although  this  technique  has 
been  successfully  applied  to  several  practical  problems  in  duct 
wave  propagation,"*  it  has  the  following  limlations:  I)  the 
series  expansions  oflen  involve  special  functi.sns  which  are 
difficult  to  compute;  2)  at  high  frequencies  and  at  the 
boundaries  the  series  are  slowly  convergent — therefore,  a 
large  number  of  terms  in  the  series  must  be  retained  to  ensure 
accurate  results,  which  often  requires  excessive  computation 
lime;  3)  this  method  can  only  be  used  with  special  coordinate 
systems  and  boundary  conditions  for  which  the  separation  of 
variables  can  be  applied.  At  present  only  eleven  suitable 
coordinate  systems  are  known  (Ref.  I , p.  5 13f0. 

Tor  arbitrarily  shaped  bodies,  the  differential  form  of  the 
Helmholtz  equation  can  be  solved  by  writing  Ihe  equation  in 
terms  of  finite  differences  (Ref.  I,  p.  703fO.  Unlike 
separation  of  variables,  this  technique  is  not  limited  to  ducts 
with  simple  geometries.  A typical  application  of  finite  dif- 
ferences is  given  by  Wynne  and  Plumblee'"  who  solved  for 
the  transverse  eigenvalues  and  eigenfunctions  of  an  annular 
duct  with  lined  walls.  This  technique  involves  the 
simultaneous  solution  of  Ihe  acoustic  potential  value  at  every 
point  within  the  duct.  Once  the  potential  values  are  known, 
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the  acoustic  p'-essure  and  velocity  can  then  be  determined.  To 
obtain  sufficient  accuracy,  fine  grid  sizes  must  be  used  which 
necessitates  large  computer  storage  requirements.  This 
drawback  was  noted  by  Baumeister, ' ' Baumeister  and  Rice,  ’’ 
and  Alfredson"  who  used  this  technique  in  studies  of  duct 
wave  propagation.  Because  of  the  storage  requirements  this 
technique  has  mainly  been  applied  to  two-dimensional 
problems.  For  three-dimensional  problems  numerical 
methods  capable  of  handling  large  matrices  must  be  used 
which  require  considerable  computer  time  and  computational 
effort.''*  This  technique  is  also  impractical  in  radiation 
problems  which  involve  infinite  domains. 

To  avoid  the  limitations  of  the  differential  formulation,  the 
integral  approach  is  employed  in  this  study.  The  integral 
approach  has  been  successfully  applied  to  a wide  range  of 
acoustic  problems.  In  determining  the  sound  radiation  field 
from  vibrating  surfaces,  integral  techniques  have  been  widely 
used.'*  '"  For  example,  Chen  and  Schweikert employed 
this  method  to  determine  the  radiation  sound  patterns  for 
three-dimensional  shapes  with  mixed  boundary  conditions. 
To  check  the  accuracy  of  the  results,  they  computed  the 
radiated  field  produced  by  a piston  vibrating  on  a sphere.  For 
this  problem  an  exact  solution  exists ' and  compares  favorably 
with  the  numerical  results.  The  integral  fonnulation  is  also 
used  to  solve  the  problem  of  scattering  by  arbitrary 
shapes.'**'  Banaugh  and  Goldsmith,  for  example,  used  this 
technique  to  investigate  the  effect  of  surface  shape''*  on 
scattered  sound  fields.  By  applying  this  method  to  a circular 
cylinder,  for  which  exact  solutions  are  available,’ '*  and 
comparing  the  exact  and  numerical  solutions,  Banaugh  and 
Goldsmith  demonstrated  the  accuracy  of  the  integral  solution 
scheme.  Although  this  method  is  capable  of  handling  mixed 
boundary  conditions,  only  surfaces  with  rigid  boundaries 
were  considered  in  Ref.  19.  The  effect  of  mixed  boundarv 
conditions  was  included  in  studies  by  Liu  and  Martenson** 
and  Quinn^’  of  the  internal  acoustic  pattern  of  lined  ducts 
with  arbitrary  shapes.  Comparison  of  the  theoretical 
predictions  with  experimental  data  showed  generally  good 
agreement.  Unpublished  work  by  Zinn  and  Gaylord *'* 
demonstrated  the  applicability  of  the  integral  formulation  for 
the  determination  of  the  natural  frequencies  and  modes  for 
two-dimensional  shapes.  In  this  study  the  accuracy  of  the 
technique  was  determined  by  comparing  the  natural 
frequencies  and  mode  shapes  with  available  exact  solutions 
for  a two-dimensional  cylinder  with  rigid  walls.  The 
agreement  is  to  within  four  decimal  places  which  is  two- 
orders-of-magnitude  more  accurate  than  previous  results 
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obtained  by  solving  the  differential  Helmholtz  equations 
using  finite  differences.'"  In  another  study  by  Tai  and 
Shaw,-'  the  integral  method  was  applied  to  a right  triangle. 
The  resulting  eigenfrequencies  compared  with  exact  solutions 
to  within  5'’7o  and  the  maximum  deviation  between  the 
numerically  computed  and  exact  potential  fields  was  less  than 
lo/o. 

To  demonstrate  the  accuracy  and  the  versatility  of  the 
integral  solution  technique,  results  are  obtained  for  several 
acoustic  problems  involving  a variety  of  geometries.  To 
obtain  a solution,  the  integral  equation  is  first  discretized  to 
form  a system  of  algebraic  equations  which  are  then  solved 
for  the  acoustic  potential  at  discrete  points  on  the  boundary. 
From  these  values  the  rest  of  the  sound  field  is  obtained. 
Methods  for  increasing  the  numerical  accuracy  by  use  of 
Gaussian  quadrature  and  other  numerical  integration 
methods  are  presented  and  discussed.  The  first  problem 
considered  is  the  numerical  evaluation  of  the  resonant 
frequencies  and  natural  modes  of  two-dimensional  circular, 
rectangular,  and  star  configurations.  Exact  and  numerical 
values  are  compared  for  the  circle  and  rectangle.  The  next 
problem  considered  is  a two-dimensional  duct  with  a right- 
angle  bend  with  a sound  source  at  one  end  and  sound  ab- 
sorption treatment  at  various  locations  along  the  duct.  The 
results  are  compared  with  finite  difference  solutions.  These 
studies  demonstrate  the  applicability  of  the  integral  for- 
mulation to  complicated  geometries  and  general  boundary 
conditions.  The  next  problem  considered  is  the  two- 
dimensional  radiation  problem  of  a piston  set  in  a right 
circular  cylinder.  Again,  the  exact  and  numerical  acoustic 
fields  are  computed  and  compared.  Finallv.  a three- 
dimensional  problem  of  determining  the  acoustic  properties 
of  a sphere  is  considered.  The  internal  field  is  obtained  using 
an  axisymmetric  formulation. 


11.  Governing  Equations 

The  integral  formulations  of  the  wave  equation  for  internal 
and  radiation  acoustic  problems  are  developed  ir.  this  section 
for  two  and  three  dimensions.  The  boundary  conditions 
generally  encountered  in  practical  problems  are  then 
discussed.  For  clarity,  only  a brief  account  of  the  derivation 
of  the  basic  equations  will  be  given  in  this  section.  For  a more 
detailed  and  rigorous  development.  Refs.  26  through  29  can 
be  consulted. 

Assume  a frictionless,  homogeneous  gas,  and  let  and  p,, 
be  the  density  and  pressure  of  the  fluid  at  rest.  Representing 
the  acoustic  pressure  and  particle  velocity  at  a time  I by  p and 
u,  Euler’s  equation  for  the  conservation  of  momentum  gives 

(1U 

p„~-  +Vp  = 0 (I) 


The  continuity  equation  yields  the  relationship 
dp  , 

^ +p„c{,'7  u^O  (2) 

0( 

where  c„  is  the  speed  of  sound.  By  defining  an  acoustic 
potential  function  ♦ such  that 

u=V+  (3t 

Equation  ( I ) provides  the  relation 

3 + 

(4) 


and  Eq.  (2)  results  in  the  classical  wave  equation 


I 


0 


(5) 


The  wave  equation  can  also  be  written  in  terms  of  p and  u,  but 
it  is  more  convenient  to  work  with  an  acoustic  potential 
function,  from  which  both  the  acoustic  pressure  and  particle 
velocity  can  readily  be  obtained. 

Equation  (5)  is  the  wave  equation  for  a general  time 
dependence  and  can  be  written  in  integral  form  and  solved  by 
using  retarded  potentials.’''^*  However,  for  most  practical 
problems  a sinu.soidal  time  dependence  can  be  assumed  which 
simplifies  the  problem  considerably.  Assume 

'Kr./)=«(r)e'“'  (6) 

Substituting  Eq.  (6)  into  Eq.  (5)  gives  the  Helmholtz  equation 

V + k^(t>  ~0  (k  = ip/Co)  (2) 

which  can  be  solved  by  simpler  methods  not  involving  the  use 
of  retarded  potentials. 


Integral  Formulation 

To  obtain  an  integral  formulation  of  the  Helmholtz 
equation,  consider  the  problem  shown  in  Fig.  1 . Applying 
Green’s  theorem  to  the  Helmholtz  equation ' gives  the 
following  integral  relation 


dCiP.O)  dibiO) 

'^(Q)  — -CiP.Q)  -^)dSy  = 0 (8) 


dn, 


dn. 


where  <t>  is  the  acoustic  potential  function  and  G is  the  Green’s 
function  defined  by  Eqs.  (14-16),  which  also  satisfies  the 
Helmholtz  equation.  The  Green’s  function  is  regular  inside 
the  surface  except  when  P-Q.  At  this  point  G is  singular.  To 
remove  this  singularity  from  the  integral  given  by  Eq.  (8), 
point  P is  surrounded  by  a small  sphere  or  circle  a of  radius  e. 
The  integral  will  now  include  a term  over  a which,  on  taking 
the  limit  as  f — 0,  gives 


<<.(P)=c|^(G(P,0) 


3<0(g) 

dnQ 


dGU^.Q) 

driQ 


Q 


(9) 


where  C is  /74  for  two  dimensions  and  1/4t  for  axisymmetric 
and  three-dimensional  shapes. 

From  Eq.  (9)  the  value  of  the  acoustic  potential  function  at 
any  point  P within  the  surface  can  be  determined  from  the 
boundary  values  of  the  potential  and  its  normal  derivative. 
Thus,  the  entire  wave  pattern  within  the  surface  can  be 
constructed.  For  arbitrarily  shaped  surfaces  for  which 
numerical  techniques  must  be  used  to  obtain  a solution,  Eq. 
(9)  requires  much  less  computer  storage  than  the  differential 
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formulation  given  by  Eq.  (7).  Using  Eq.  (9)  the  value  of  the 
potential  at  each  interior  point  can  be  obtained  independently, 
whereas  the  method  of  finite  differences  used  to  solve  Eq.  (7) 
requires  the  simultaneous  solution  of  <t>  for  every  interior 
point.  The  integral  formulation  avoids  the  large  matrices 
involved  w ith  finite  differences. 

If  the  values  of  both  <t>  and  h<i>lbn  are  known  at  every  point 
on  the  boundary  then  the  wave  pattern  can  readily  be 
determined  from  Eq.  (9).  However,  for  most  practical 
acoustic  problems  either  d<t>/dn  or  an  admittance  condition 
relating  <>  and  d<t>/dn  are  given.  Therefore,  the  values  of  the 
acoustic  potential  at  the  boundary  must  first  be  determined. 
The  necessary  relation  is  obtained  by  letting  the  point  t 
approach  the  boundary  at  some  point  T to  obtain  the 
following  relation’'’ 


'i>(T)=2C 


d<t>(Q)  dG(T,0) 

[G(T.Q)-^’--<t>(Q)  --^^IdSg  (10) 

ant 


dn. 


where  6 is  the  Dirac  delta  function.  The  Green’s  functions 


G[P,Q)  =Ho  {kr)  for  two  dimensions. 


(14) 


f ' e'*' 

G{P,Q)=2\ cos  m8d9  for  axisymmetric  bodies  (15) 

Jo  r 


and 

G(P,Q)  =e~*' / r for  three  dimensions  (16) 

where  r is  the  distance  between  points  P and  Q,  and 
ny'  (kr)  is  the  zeroth  order  Hankel  function  of  the  first  kind. 

Boundary  Conditions 

The  two  most  common  boundary  conditions  in  practical 
acoustic  problems  are  the  Neumann  and  Robin  conditions. 
The  Neumann  condition  of  interest  in  the  present  study  is 


Eq.  (10)  is  applicable  to  a smooth  boundary,  but  has  been 
extended  to  include  cusps  and  corners.  To  obtain  the 
interior  wave  pattern,  Eq.  (10)  is  first  solved  for  the  boundary 
values  of  <t>.  These  values  are  then  substituted  into  Eq.  (9)  to 
determine  the  acoustic  potential  at  the  interior  points.  Both 
Eqs.  (9)  and  (10)  involve  singular  integrands  as  T approaches 
Q although  for  smooth  surfaces  the  integrals  themselves  are 
regular. 

For  exterior  problems,  analogous  expressions  to  Eqs.  (9) 
and  ( 10)  are  obtained  by  taking  the  point  P outside  the  surface 
r.'*  The  integration  in  Eq.  (8)  is  then  carried  out  over  the 
boundary,  around  a circle  or  sphere  of  radius  e with  point  P as 
a center,  and  then  around  a circle  or  sphere  of  radius  /?,  which 
is  arbitrarily  large.  In  this  manner  the  integration  includes  the 
entire  external  domain.  However,  by  applying  Sommerfeld’s 
radiation  condition,  it  can  be  shown  that  the  integral  about 
the  infinite  sphere  or  circle  approaches  zero  as  P approaches 
infinity.^  Thus,  the  corresponding  equations  for  the  external 
domain  become 


f d<t>(Q)  dG{P,Q) 

4>(P)  = -C\  lG{P.Q)-^—-<t>{P)  : ^ |d5t,  (11) 

J r on,,  dn„ 


d<t>/dn=A  (17) 

where  A is  the  velocity  amplitude  of  a given  sound  source.  In 
the  absence  of  a sound  source  A=0;  this  condition  means  that 
the  particle  velocity  is  zero  at  the  boundary  which  implies  a 
perfectly  reflecting,  or  rigid  surface.  For  surfaces  which 
absorb  sound,  such  as  lined  duct  walls,'  an  admittance  con- 
dition is  usually  specified,  which  leads  to  the  Robin  condition. 
Defined  as  the  radio  of  the  normal  component  of  the  particle 
velocity  to  the  pressure  perturbation,  the  admittance  y can  be 
written  as 

y^PoCoiuJp)  (18) 

Substituting  for  u„  andp  from  Eqs.  (3)  and  (4)  gives 

(d<f>/dn)  +iky<t>  = 0 (19) 

Eq.  (19)  is  the  Robin  condition.^’  For  sound-absorbing 
materials  or  devices,  the  admittance  can  be  either  analytically 
determined or  measured  using  the  impedance  tube  or  a 
related  technique.’’  ’^  The  effects  of  a given  material  on  the 
internal  acoustic  properties  of  a particular  geometry  can  be 


and 

<t>(T)  = 


[GiT.Q) 


dn„ 


dG(T,Q) 

-0(0)----  - IdS^, 

orin 


(12) 


It  is  important  to  note  that  Eqs.  (11)  and  (12)  involve  in- 
tegrations about  the  boundary  of  the  body  only.  Thus,  the 
radiated  field  at  any  distance  from  the  body  can  be  obtained 
once  the  surface  acoustic  potential  is  known.  With  finite 
differences,  the  values  of  the  potential  at  every  point  in  a very 
large  domain  would  have  to  be  computed  in  order  to  obtain 
the  radiated  field.  Also,  an  artificial  boundary  condition  at  a 
large  distance  from  the  surface  must  be  assumed.  These 
factors  make  the  application  of  finite  differences  to  problems 
of  this  type  rather  inefficient  whereas  the  integral  formulation 
can  readily  be  adopted  to  such  situations. 

Eqs.  (9)  through  (12)  are  applicable  to  two-dimensional, 
axisymmetric,  and  three-dimensional  acoustic  problems.  In 
the  two-dimensional  and  axisymmetric  cases,  these  equations 
involve  line  integrals;  and  in  the  three-dimensional  case,  the 
integrals  are  taken  over  a surface.  Note  that  the  dimen- 
sionality of  the  problem  is  reduced  by  one — a valuable 
simplification. 

The  Green’s  functions  satisfy  the  following  inhomogeneous 
forms  of  the  Helmholtz  red'uction  with  homogeneous 
boundary  conditions ' 

v’G  + k’G  = b{P-Q)  (13) 
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determined  by  substituting  the  admittance  of  the  material  into 
Eq.  (19)  and  solving  Eqs.  (9)  and  (10)  or  (1 1)  and  (12)  for  the 
acoustic  potential.  Thus,  the  analytical  techniques  used  in  this 
investigation  is  applicable  to  a vast  number  of  duct  acoustic 
problems.  Since  the  admittance  of  a combustion  process  can 
also  be  measured,  this  analysis  can  also  be  applied  to  related 
linear  combustion  instability  problems,  provided  that  the 
equations  are  applied  to  regions  where  the  Helmholtz 
equation  holds  and  mean  flow  effects  can  be  neglected.  By 
replacing  the  combustion  process  by  an  admittance  condition, 
studies  of  combustion  instability  have  been  conducted  in 
liquid  and  solid  propellant  combustors.  This  research 
allows  the  extension  of  these  analyses  to  more  general  shapes. 

Substituting  Eq.  (19)  into  Eq.  (10)  gives  for  the  internal 
field 

f dC(T,Q) 

(i.(D+2C  »(Q)[  ' — +iky(Q)G{T.Q)]dSQ 

J r driQ 

= - 2cj  ^ A {Q)G(  T. Q) dSg  (20) 

\ similar  expression  is  obtained  from  Eq.  (12)  for  the  exterior 
problem.  For  surfaces  with  spatially  varying  admittances,  the 
admittance  is  a function  of  Q.  For  most  cases  considered  in 
this  study,  y is  assumed  constant  although  nonuniform  ad- 
mittance distributions  can  be  easily  handled. 

III.  Solution  Technique 

In  the  last  section,  the  integral  equations  were  developed 
which  describe  the  interior  and  exterior  acoustic  fields  or  a 
surface  with  arbitrary  shape  and  mixed  boundary  conditions. 
The  numerical  solution  technique  for  solving  these  equations 
to  obtain  the  internal  or  radiated  acoustic  patterns  is 
presented  in  this  section  and  can  be  divided  into  four  parts. 
The  first  is  the  discretization  of  the  integral  equations  into  a 
corresponding  system  of  linear,  and  algebraic  equations  in 
suitable  for  solution  on  a computer.  The  second  part  is  the 
specification  of  the  geometry  and  boundary  conditions.  The 
third  is  the  computation  of  the  coefficients  of  the  system  of 
equations  and  the  final  part  is  the  methods  used  to  solve  for 
the  surface  potential  from  the  algebraic  equations. 

Diserelizalion  of  the  Integral  Equations 

In  two  dimensions  and  for  axisymmetric  problems,  Eq.  (20) 
involves  1 one-dimensional  improper  integral  about  the 
boundary  line.  For  this  type  of  problem  several  numerical 
integration  techniques'"  '^  are  available.  The  simplest  is  the 
trapezoidal  rule  which  has  been  shown  to  yield  excellent 
results  in  two-dimensional  studies  with  this  type  of  in- 
tegral.Using  this  numerical  integration  scheme,  Eq. 
(20)  becomes 

<t>„.y2cYd  + ikyiG(.r„„)\^S, 

, , 6n,  J 


= -2cl^A,G{r„„)^S,  (21) 

, I 

A here  one  equation  for  (t>  is  obtained  for  each  value  of  m and 
m is  varied  from  1 to  N.  Eq.  (21)  was  initially  used  in  this 
investigation  to  generate  the  N equations  for  <t>  and  accurate 
results  were  obtained  when  the  admittance  y was  zero 
everywhere  on  the  boundary"”  which  is  the  case  considered  in 
previous  studies.  However,  when  a nonzero  admittance 
IS  assumed,  this  technique  gives  inaccurate  results  because  of 
the  contribution  from  the  Green’s  function  when  the  pointy 
appr  taches  m.  Because  of  the  singular  nature  of  the  Green’s 
function  at  the  point  m,  care  must  be  taken  when  numerically 
integrating  this  function  over  the  subinterval  m.  To  increase 
the  accuiacy  in  evaluating  the  integrand,  Eq.  (20)  is  broken  up 
into  N integrals  given  by  Eq.  (22). 


4>„ 


+ </i 

■Sm  - 


dG(r,^) 

dzi'i 


+ iky, 


,C(r,^)j  dsj 


+ 2C 


S' 


E 

J=l 


j*m 


G(rj„)dSj 


(22) 


In  both  Eqs.  (21)  and  (22)  the  values  of  <t>  are  assumed  to  be 
constant  over  each  of  the  N subintervals.  The  difference  is  the 
method  by  which  the  terms  involving  the  Green’s  function  are 
evaluated.  In  Eq.  (21)  an  average  value  is  computed  over  each 
of  the  subintcrvals  based  on  With  Eq.  (22)  these  terms  are 
integrated  numerically  from  „ to  G + 'j.m  using  Gaussian 
quadrature'”'"®  to  obtain  more  accurate  values.  This  type  of 
formulation  has  been  used  before  with  trapezoidal  instead  of 
Gaussian  quadrature  formulas.'"  In  the  present  study  for 
two-dimensional  and  axisymmetric  problems,  a reduction  in 
error  of  two  orders  of  magnitude  in  the  numerical  results  for  a 
nonzero  admittance  was  achieved  using  Eq.  (22)  instead  of 
Eq.(2l)."” 


Surface  (iromciry  and  Boundary  Conditions 

The  firs!  step  in  solving  Eq.  (22)  is  the  determination  of  the 
coefficients  of  0 and  <t>„ . These  coefficients  depend  upon  the 
surface  geometry  through  the  terms  d/drij,  and  AS^.  By 
specifying  the  admittance  y and/or  the  sound  velocity  am- 
plitude A over  every  subinterval  j,  the  effect  of  the  boundary 
conditions  are  included  in  the  evaluation  of  the  coefficients. 

To  solve  for  the  terms  involving  the  surface  geometry,  the 
first  expression  inside  the  integrals  of  Eq.  (22)  is  written  as 

dG{r)  _ dG(r)  dr 
dn  dr  dn 

Th  expressions  for  dG/dr  are  obtained  by  differentiating  Eqs. 
(14)  through  (16).  Substituting  this  expression  into  Eq.  (22) 


Table  I Ggenfrequencies  and  natural  modes  of  a circle 
for  various  admittance  values 


ADMITTANCE  VALUE 


■'•ov 


o 


\ COMPUTED  1 

1 64122 

183244-0  4423^ 

14441-0  0071 

1 EXACT  \ 

184118  ! 

l8322  + 0.443a 

14384 

. COMPUTED [ 

9 05423 

sOTaiTOSsSTi! 

25369'00f5ll 

/ EXACT  i 

3 09424 

|30796«0  5A«( 

2 5427 

.computed! 

3 83175 

jNot  Coinputtd 

, Not  Computod 

/ EXACT  ' 

3 83I7I 

',3,8188^0  30961 

3 9510 

\ COMPUTED! 

420135 

42538  061991 

39816 -0  0231 

'exact 

4 20119 

4 253  2 0 6301, 

39619 

\ COMPUTED 

93IT83  1 

InoI  Cofnputtd 

Not  Computod 

/ EXACT 

531759 

^3953  0 71011 

45767 

Table  2 Resonant  f requencies  and  natural  modes  of  a 
rectangle  for  different  admittance  values  at  the  ends 


ACMITTANCE  VALUE 


r-o 

y • 03 

y - OSl 

3 1432 

3 1416 

3190  4-061991 

3 142  4 0 61901 

2 55840  002 1 

2 999 

6 2877 

6 2832 

6 30240  61961 

{6  283  40  61901 

5 686  40  0021 

98i» 

7 0312 

7 0246  1 

7 146  40988 1| 
Not  Cotnputod 

6 333  - 0 0071 

6 283 

8 893 

6 886 

8 934  40  6101  | 
Not  Cofnputod 

8 303*0011) 

8299 

9 4329  j 

9 4248 

9 456  4 061061 

9 429  4 0 6M01 

8 047  4 0 007t 

6*421 

uttq 

EXACT  j 
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gives  a relation  which  involves  dridnj,  rj„  and  dS.  The  ex- 
pressions for  dr/drij,  and  dS  can  be  written  in  parametric 
form  for  bodies  with  simple  shapes,  and  this  type  of 
representation  has  been  used  in  previous  studies  using  simple 
geometries. By  taking  advantage  of  symmetry, 
considerable  savings  in  computer  storage  and  computation 
times  were  achieved.  In  fact,  Greenspan  and  Werner^’ 
showed  that  for  a circle,  Eq.  (22)  can  be  reduced  to  a single 
equation  instead  of  a system  of  equations  which  could  readily 
be  solved  to  obtain  the  acoustic  field.  In  the  study  by  Tai  and 
Shaw,^’  the  method  of  images  (Ref.  1,  Ch.  11)  was  used  to 
greatly  reduce  the  numer  of  points  necessary  to  compute 
eigenfrequencies  and  eigenmodes  of  a family  of  triangles. 
Although  these  studies  demonstrate  valuable  simplifications 
which  can  be  made  in  applying  the  integral  formulation  to  a 
particular  problem,  the  techniques  used  are  not  applicable  to 
more  general  problems  involving  complicated  geometries  and 
nonuniform  boundary  conditions. 

In  the  present  study  for  two  dimensions,  the  expressions  for 
the  geometric  variables  are  written  in  parametric  form  only 
for  the  circle.  In  the  rest  of  the  configurations  considered,  a 
general  formulation  is  used.  The  fact  that  a parametric 
representation  cannot  be  used  in  general  cases  is  not  a serious 
drawback— in  fact,  it  somewhat  simplifies  the  formulation. 
Consider  the  general  two-dimensional  problem  depicted  in 
Fig.  2.  By  specifying  the  x and  y coordinates  at  the  midpoint 
of  each  of  the  subintervals,  the  distance  is  readily  com- 
puted from  the  expression 

rjm  = '^(Xj-x„)^+ (yj-yJ^  (23) 

The  expression  for  dr/duj  can  then  be  obtained  since  it 
represents  the  dot  product  of  the  gradient  of  r and  the  normal 
at  j.  Thus, 

_ (x,-x„)n,,+  (yj~y„)ny, 

dn,  rj„ 

where  is  the  component  of  the  normal  vector  j in  the  x 
direction  (or  the  cosine  of  the  angle  between  the  normal  vector 
and  the  x-axis)  and  riyj  is  the  corresponding  > component  (the 
sine  of  the  angle  between  the  normal  vector  and  the  j'-axis). 
Analagous  expressions  for  r,„  and  dr/drij  can  be  obtained  for 
axisymmetric'’  '*  and  three-dimensional  problems.  For  two- 
dimensional  and  axisymmetric  problems,  the  line  segment 
length  ^S,  is  simply 

ixSj  = si  {Xj  * -Xj.„,)^+  {yj^,,-y,-,„)^ 

or,  for  (V equally  spaced  subintervals,  Sj  = I /N  where  L is  the 
length  of  the  perimeter  of  the  surface.  For  three-dimensional 
bodies,  AS,  is  the  area  of  each  of  the  subsurfaces  taken  over 
the  boundary. 


Compulation  of  the  Coefficients  of  the  Discretized  Integral  Equation 

Once  the  geometry  has  been  specified,  the  coefficient  of  <t> 
in  Eq.  (22)  can  be  determined  by  evaluating  the  Green’s 
functions  C(/-j„)  and  3G(/-^„) There  are  two  problems 
in  determining  these  functions:  the  first  is  a rapid,  accurate 
method  for  computing  them  over  a wide  range  of  the 
argument  krj„;  and  the  second  is  the  singularity  associated 
with  each  function  as  approaches  zero. 

For  the  two-dimensional  problems  to  compute  the  Hankel 
functions  two  routines  have  been  used  in  this  study.  The  first 
consists  of  a series  expansion  using  standard  formulas  for  the 
Hankel  function  with  complex  arguments.  A sufficient 
number  of  terms  is  taken  to  satisfy  a specified  degree  of 
accuracy.  To  minimize  time,  a different  series  expansion 
which  was  developed  by  Hitchcock”  is  used  for  determining 
these  functions  in  the  studies  of  the  rectangle,  star,  and  duct 
with  a right-angle  bend.  With  his  formulation,  accuracies  of 
10  or  greater  are  achieved  using  nine  terms  or  less  in  the 
series  expansion.  Reductions  of  up  to  50%  in  computer  times 
can  be  achieved  with  this  formulation. 

For  the  axisymmetric  problem,  the  integral  in  Eq.  (15)  is 
carried  out  using  a 20-point  Gauss- Legendre  quadrature 
formula.  For  three-dimensional  problems,  evaluation  of  the 
Green’s  function  given  by  Eq.  (16)  is  straightforward. 

The  major  problem  in  accurately  computing  the  coef- 
ficients in  Eq.  (22)  is  the  singularity  associated  with  the 
Green’s  functions  as  rj„  approaches  zero;  that  Is,  as  the  point 
j approaches  m in  Fig.  2.  The  two-dimensional  and 
axisymmetric  Green’s  functions  have  logarithmic 
singularities.  In  this  study,  the  inaccuracies  involved  are 
minimized  by  subdividing  the  intervals  as  indicated  by  Eq. 
(22). 

Determination  of  the  Acoustic  Potential 

Once  the  coefficients  of  the  surface  potential  at  each 
discrete  point  on  the  surface  are  determined,  the  equations  are 
solved  for  <j>  using  a complex  Gauss-Jordan  reduction  scheme. 
The  interior  or  exterior  points  can  then  be  found  using  the 
discretized  form  of  Eq.  (9). 

To  determine  the  eigenfrequencies  of  a particular  geometry, 
the  technique  described  in  Ref.  40  is  used.  Essentially,  this 
technique  consists  of:  .1)  determining  the  frequency  k for 
which  the  determinant  of  the  coefficients  in  the  homogeneous 
form  of  Eq.  (22)  is  zero,  2)  normalizing  the  equation  at  the 
eigenfrequency  to  obtain  the  surface  distribution  of  the  mode, 
and  3)  using  Eq.  (9)  in  discretized  form  to  find  the  interior 
sound  field. 


V.  Results 

Using  the  numerical  techniques  described  in  the  last  section, 
solutions  have  been  obtained  for  a variety  of  two-dimensional 
and  axisymmetric  problems  to  demonstrate  its  broad  range  of 
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Fig.  3 Nodal  points  and  lines  for  the  first  nine  modes  of  the  star. 


Table  3 Surface  potentials  for  a circle  of  unit  radius. 


second  mode,  1=0 

Angle 

Numerical 

Evaci 

12 

.91.30 

24 

.6681 

.6691 

36 

..3077 

.3090 

48 

- .1059 

- .I04S 

60 

- .5012 

50(X1 

72 

- .8099 

.8090 

84 

9785 

- .9181 

96 

- .9779 

- .9’8I 

108 

- .8082 

.8090 

120 

- .4988 

.5(XX) 

132 

- .1013 

- .1045 

144 

.3104 

.3090 

156 

.6702 

.6691 

168 

.9141 

.9135 

180 

1.0000 

i.onoo 
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Table  4 Surface  potential  fora  rertanxie,  height  to  width 
ratio  = O.S,  first  mode,  rigid  walls 


y 

Numerical 

Exact 

1/2 

1/14 

1.0061 

1.0063 

12 

2/14 

1.0060 

1.0063 

1/2 

.1/14 

1.0055 

1.0063 

1.1/28 

1/4 

1.0000 

1.0000 

1 1/28 

1/4 

.9501 

.9499 

9/28 

1/4 

.8524 

.8521 

7 '28 

1/4 

.7119 

.7116 

5/28 

1/4 

.5156 

.5354 

.V28 

1/4 

.1325 

.3324 

1/28 

1/4 

.1127 

.1127 

Fig.  S Comparison  of  numerical  results  for  a duct  with  a right-angle 
bend  using  the  integral  and  finite  difference  approaches.  Case  I. 


applications.  The  two-dimensional  form  of  the  integral 
equation  has  been  used  to  compute  the  -csonant  frequencies 
and  natural  modes  of  a circle,  rectangle,  and  star  con- 
figuration. In  addition,  the  problem  of  a duct  with  a right 
angle  bend  is  considered,  and  results  using  Eq.  (22)  are 
compared  with  finite  difference  solutions.  The  two- 
dimensional  problem  of  sound  radiation  from  a right  circular 
cylinder  is  then  considered  and  the  numerical  and  exact 
solutions  are  compared.  Finally  the  acoustic  properties  of  a 
sphere  are  computed  using  the  axisymmetric  formualtion. 

For  a circle  and  rectangle,  comparisons  between  exact  and 
numerical  solutions  are  presented  in  Tables  1 and  2.  In  these 
tables  the  numerical  and  exact  eigenfrequencies  are  tabulated 
for  three  admittance  values,  >"=0,  y = 0.3,  >>  = 0.3/,  with  thirty 
points  taken  on  the  boundary.  The  best  agreement  between 
the  computed  and  exact  results  occurs  at  the  zero  admittance 
condition.  For  the  circle,  the  real  part  of  the  eigenfrequencies 
compare  to  five  significant  figures  and  the  imaginary  parts  are 
accurate  to  0.001  for  the  first  five  modes.  When  a nonzero 
admittance  condition  is  introduced,  the  accuracy  is  reduced  to 
three  significant  figures  in  the  real  part  and  to  0.01  in  the 
imaginary  part  of  the  eigenfrequencies. 

As  with  the  circle,  the  agreement  between  the  exact  and 
numerical  values  for  the  rectangle  is  good  for  a rigid  bound- 
ary but  deteriorates  when  a nonzero  admittance  is  in- 
troduced. From  Table  2 the  agreement  is  to  almost  four 
significant  figures  in  the  real  part  of  the  eigenfrequency  and  to 
within  0.01  in  the  imaginary  part  for  a rigid  wall.  The 
Gaussian  integration  techniques  developed  in  Sec.  Ill  improve 


y.| 

Fig.  4 I.«calian9  of  the  discrete  points,  nonzero  admittim 
boundnries,  and  the  sound  source  for  the  duct  with  a right-angie  ben* 


Fig.  6 Comparison  of  numerical  results  for  a duct  with  a right-angle 
bend  u.sing  the  integral  and  finite  difference  approaches.  Case  2. 


the  accuracy  of  the  computed  eigenfrequencies  for  a nonzero 
admittance  condition  by  an  order  of  magnitude. 

For  the  circle  the  accuracy  of  the  computed  natural  mode 
shapes  is  shown  in  Table  3.  The  agreement  between  the  exact 
and  computed  eigenmodes  for  a rigid  boundary  is  to  within 
0.01  “7o  for  interior  points  sufficiently  far  removed  from  the 
boundary.  For  a nonzero  admittance  at  the  surface,  the  ac- 
curacy is  to  within  2%.  These  results  are  obtained  using  the 
interior  analog  of  Eq.  (21)  which  explains  the  deterioration  in 
accuracy  of  the  interior  points  as  the  boundary  is  approached. 
Equation  (22)  is  used  in  the  studies  of  the  rectangle,  star,  and 
duct  problems  and  more  accurate  results  are  obtained  close  to 
the  boundary.  For  the  rectangle,  the  boundary  values  of  the 
acoustic  potential  are  presented  in  Table  4.  The  agreement 
between  the  exact  and  numerical  results  is  within  one-half  of  a 
percent.  Computation  times  range  from  ten  sec  per  eigen- 
frequency for  the  circle  to  45  sec  for  the  rectangle  on 
UNIVAC  1 108  computer.  Using  the  discretized  form  of  Eq. 
(9),  interior  points  require  approximately  two  sec  per  point  to 
compute. 

In  studying  the  star-shaped  boundary,  which  is  of  interest 
in  solid-rocket  combustion  instability  problems,  the  ap- 
plicability of  the  integral  solution  technique  to  a complicated 
geometry  for  which  separation  of  variables  does  not  apply  can 
be  assessed.  The  first  nine  eigenfrequencies  and  natural  modes 
for  the  star  are  presented  in  Fig.  3 for  a rigid  wall  with  48 
points  taken  on  the  surface.  The  most  unique  feature  of  the 
acoustic  field  for  the  star  is  the  appearance  of  nodal  points  at 
some  of  the  resonant  modes.  In  the  circle  and  rectangle  nodal 
lines  only  are  present,  and  they  follow  one  of  the  separable 
coordinates  of  the  boundary.  With  the  star  both  nodal  lines 
and  points  can  occur  which  is  in  qualitative  agreement  with 
experimental  observations  for  unstable  solid  propellant 
combustors.  Computation  times  are  from  60  to  73  sec  per 
mode.  The  modes  of  a typical  solid  propellant  configuration 
during  a burn  have  also  been  computed  and  are  given  in  Ref 
45. 

The  last  internal  two-dimensional  problem  investigated  is 
that  of  a duct  with  a right-angle  bend  shown  in  Fig.  4.  The 
reasons  for  studying  this  configurationarejJ)^to_investig^(Sj. 
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nonuniform  surface  admittance,  2)  to  include  a sound  source 
in  the  integral  formulation,  and  3)  compare  the  results  ob- 
tained by  the  integral  technique  with  the  finite  difference 
solutions  of  Ref.  13. 

The  results  obtained  using  this  configuration  are  presented 
in  Figs.  5 and  6 and  are  compared  with  the  solutions  obtained 
uisng  the  finite  difference  method.  Although  the  results  using 
the  integral  approach  are  in  qualitative  agreement  with  the 
finite  difference  solution,  quantitative  agreement  is  lacking. 
The  same  number  of  boundary  points  are  taken  in  both  cases. 
Doubling  the  number  of  subintervals  using  Eq.  (22)  does  not 
improve  the  agreement  between  the  two  sets  of  data. 
However,  it  does  show  that  the  results  of  the  integral  for- 
mulation are  self-consistent.  An  experimental  setup  is 


I able  S Kesiinani  f requencies  and  surface  potentials  for  a sphere 
of  unit  radius,  first  and  second  modes,  axisymmetric  formulation 


( oinpiiicd 

Resonant  frequencies 

Exact 

2.084-0.(K)41 

2,082 

3..346-0.(K)7i 

3.342 

Normalized  surface  potential 

Angle 

Computed 

Exact 

(First  mode) 

s 

1.00000 

1.00000 

15 

96960 

.96%2 

25 

.90975 

.90977 

35 

.82232 

.82228 

45 

.70992 

.70981 

ss 

.57593 

.57577 

65 

.42439 

.42423 

75 

.2.5992 

.25981 

85 

.08753 

.08749 

(Second  mode) 

S 

1.0000 

1.0000 

15 

,9(»K 

9099 

25 

.7403 

.7405 

35 

.5121 

.5124 

45 

.2275 

.2529 

55 

-.0072 

- .0066 

65 

- .2346 

- .2348 

75 

-.4051 

- .4041 

85 

- .4953 

-.4942 

‘Numerkal  re^ullsnbiamed  using  Fq.  (21)  instead  of  Fq.  (22|. 


currently  being  developed  to  check  these  results  and  should 
clarify  the  discrepancy  between  these  two  methods. 

For  two-dimensional  radiation  problems,  excellent  results 
are  obtained  as  shown  in  Fig.  7.  Here  the  radiated  field  from  a 
piston  set  in  a right  circular  cylinder  is  computed  and  com- 
pared with  exact  results  from  Ref.  3.  The  mean  square  error  is 
less  than  2%  while  the  computation  time  required  is  15  sec  to 
obtain  both  the  surface  and  far  field  patterns. 

To  check  the  accuracy  of  the  axisymmetric  formulation,  the 
first  two  resonant  frequencies  and  natural  mofie  shapes  of  a 
sphere  were  computed  and  are  presented  in  Table  5.  As  with 
the  two-dimensional  problems,  agreement  between  the  exact 
and  numerical  calculations  is  excellent.  Computation  times 
are  approximately  two  minutes  per  mode;  however,  no  at- 
tempt was  made  to  take  advantage  of  the  symmetry  of  the 
problem  which  can  reduce  the  computation  time  by  at  least  a 
factor  of  two. 

Conclusions  and  Recommendations 

The  results  for  the  circle  and  rectangle  show  that  the  in- 
tegral technique  is  very  accurate  in  determining  resonant 
frequencies  and  natural  mode  shapes.  Its  application  to  the 
star  configuration  demonstrates  its  usefulness  in  studying  the 
acoustics  of  complicated  shapes.  For  the  duct  with  a right- 
angle  bend,  the  integral  approach  is  shown  to  be  applicable  to 
nonuniform  boundary  conditions  involving  sound  sources. 
The  formulation  also  gives  accurate  results  for  two- 
dimensional  radiation  problems  shown  in  the  study  of  the 
right  circular  cylinder. 

With  the  axisymmetric  formulation  accurate  results  are 
obtained  for  the  internal  eigenmodes  of  a sphere.  Extensions 
to  more  complicated  boundaries  can  readily  be  made. 
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■ 

This  report  summarizes  the  work  performed  during  the  first  year  of  a research  effort  to 
determine  the  sound  fields  associated  with  jet  engine  inlet  configurations.  A solution 
approach  for  axisymmetric  bodies  based  upon  the  integral  formulation  of  the  wave  equation 
has  been  developed.  This  solution  approach  circumvents  the  uniqueness  problems  which 
normally  occur  at  certain  frequencies  when  "straight  forward"  solutions  of  the  integral 
equation  are  obtained.  A numerical  method  and  a computer  program  for  solving  for  the 
acoustic  field  associated  with  general  inlet  confi'gurations  and  boundary  conditions  have 
also  been  developed.  To  evaluate  the  numerical  method,  computed  and  exact  results  ares^ 
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compared  for  a sphere  and  a finite  length  cylinder.  For  continuous  boundary  conditions, 
the  agreement  is  within  ten  per  cent  over  a range  of  nondimensional  frequencies  from  one 
to  ten.  For  discontinuous  boundary  conditions,  the  numerical  errors  increase  by  a 
factor  of  two.  This  report  presents  results  for  a given  inlet  configuration  and  the 
computed  and  exact  solutions  are  shown  to  agree  to  within  ten  per  cent  over  the  non- 


dimensional  frequency  range  from  one  to  ten. 
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